Chao
The chao calculator returns the Chao1 richness estimate for an OTU definition. This calculator can be used in the summary.single, collect.single, and rarefaction.single commands. The calculations for the bias-corrected Chao1 richness estimator are implemented as described by Chao in the user manual for her program spade and modified by Colwell in his program estimates.
\[S_{chao1} = S_{obs} + \frac{ {n_1}\left ({n_1}-1 \right )}{2\left ({n_2}+1 \right )}\]where,
\(S_{chao1}\) = the estimated richness
\(S_{obs}\) = the observed number of species
\(n_{1}\) = the number of OTUs with only one sequence (i.e. “singletons”)
\(n_{2}\) = the number of OTUs with only two sequences (i.e. “doubletons”)
To calculate the 95% confidence intervals we assume a lognormal distribution of the variance:
\(var\left ( S_{chao1} \right ) = \frac{ {n_1}\left ({n_1}-1 \right )}{2\left ({n_2}+1 \right )} + \frac{ {n_1}\left (2{n_1}-1 \right )^2}{4\left ({n_2}+1 \right )^2} + \frac{ {n_1}^2{n_2}\left ({n_1}-1 \right )^2}{4\left ({n_2}+1 \right )^4}\), when n1>0 and n2>0
\(var\left ( S_{chao1} \right ) = \frac{ {n_1}\left ({n_1}-1 \right )}{2} + \frac{ {n_1}\left (2{n_1}-1 \right )^2}{4} - \frac{ {n_1}^4}{4S_{chao1}}\), when n1>0 and n2=0
\(var\left ( S_{chao1} \right ) = S_{obs} e^{\left (-N / S_{obs} \right )}\left (1- e^{\left (-N / S_{obs} \right )}\right )\), when n1=0 and n2>=0
If n1>0, then:
\[C = exp \left (1.96 \sqrt{\ln \left ( 1 + \frac{var\left ( S_{chao1} \right )}{\left ( S_{chao1} - S_{obs}\right )^2 }\right )} \right )\] \[LCI_{95\%} = S_{obs} + \frac {S_{chao1} - S_{obs}}{C}\] \[UCI_{95\%} = S_{obs} + C \left ( {S_{chao1} - S_{obs}} \right )\]Otherwise:
\[P = e^{\left (-N/S_{obs}\right)}\] \[LCI_{95\%} = max \left (S_{obs}, \frac{S_{obs}}{1-P}-1.96 \left ( \frac{S_{obs}P}{1-P} \right)^\frac{1}{2}\right )\] \[UCI_{95\%} = \frac{S_{obs}}{1-P} + 1.96 \left ( \frac{S_{obs}P}{1-P} \right)^\frac{1}{2}\]where,
LCI = Lower bound of confidence interval
UCI = Upper bound of confidence interval
Open the file 98_lt_phylip_amazon.fn.sabund generated using the Amazonian dataset with the following commands:
mothur > cluster(phylip=98_lt_phylip_amazon.dist, cutoff=0.10)
The 98_lt_phylip_amazon.fn.sabund file is also outputted to the terminal window when the cluster() command is executed:
unique 2 94 2
0.00 2 92 3
0.01 2 88 5
0.02 4 84 2 2 1
0.03 4 75 6 1 2
0.04 4 69 9 1 2
0.05 4 55 13 3 2
0.06 4 48 14 2 4
0.07 4 44 16 2 4
0.08 7 35 17 3 2 1 0 1
0.09 7 35 14 3 3 0 0 2
0.10 7 34 13 3 2 0 0 3
The first column is the label for the OTU definition and the second column is an integer indicating the number of sequences in the dominant OTU. The Chao1 estimator is then calculated using the values found in the third (the number of singletons)and fourth (the number of doubletons) columns. For example, chao for an OTU definition of 0.10 would be calculated as:
\[S_{chao1} = 55 + \frac{34\left (33 \right )}{2\left (14 \right )} = 95.07\] \[var\left ( S_{chao1} \right ) = \frac{34\left (33 \right )^2}{2 \left(14 \right)} + \frac{34\left (2 \left(34 \right)-1 \right )^2}{4\left (14 \right )^2} + \frac{34^2\left ( 13 \right ) \left (33 \right )^2}{ 4\left (14 \right )^4} = 341.2494\] \[C = exp \left (1.96 \sqrt{\ln \left ( 1 + \frac{341.2494}{\left (95.07 - 55\right )^2 }\right )} \right ) = 2.3641655\] \[LCI_{95\%} = 55 + \frac {95.07 - 55}{2.3641655} = 71.9\] \[UCI_{95\%} = 55 + 2.3641655\left ( {95.07 - 55} \right ) = 149.7\]Running...
mothur > summary.single(calc=chao)
...and opening 98_lt_phylip_amazon.fn.summary gives:
label Chao Chao_lci Chao_hci
unique 1553.000000 658.490667 3870.016393
0.00 1141.500000 522.185603 2658.668444
0.01 731.000000 376.708692 1527.725165
0.02 1251.000000 533.812860 3124.532743
0.03 480.428571 262.946610 962.226261
0.04 315.600000 192.960027 572.578658
0.05 179.071429 123.996099 293.627620
0.06 143.200000 103.292208 228.234805
0.07 121.647059 91.793695 186.052404
0.08 92.055556 73.303947 135.389388
0.09 96.666667 74.012067 149.489903
0.10 95.071429 71.949979 149.732827 <---
These are the same values that we found above for a cutoff of 0.10.