X-Git-Url: https://git.donarmstrong.com/?a=blobdiff_plain;f=geom.cpp;h=31580db15b7bbb08b648cd10ff3c0bcf4cd98a10;hb=e4827e0945cbda536064e5a345996b2a7dfcbb81;hp=8b484de49bcd27c772ead83f10cbf77aa14adf91;hpb=eb1c88346fb246e95a6b38935b103f95e38b82ca;p=mothur.git diff --git a/geom.cpp b/geom.cpp index 8b484de..31580db 100644 --- a/geom.cpp +++ b/geom.cpp @@ -3,12 +3,11 @@ * Mothur * * Created by Thomas Ryabin on 2/23/09. - * Copyright 2009 __MyCompanyName__. All rights reserved. + * Copyright 2009 Schloss Lab UMASS Amherst. All rights reserved. * */ #include "geom.h" -#include "calculator.h" /***********************************************************************/ @@ -43,7 +42,7 @@ RAbundVector Geom::getRAbundVector(SAbundVector* rank){ /***********************************************************************************/ EstOutput Geom::getValues(SAbundVector* rank){ try { - data.resize(2,0); + data.resize(3,0); rdata = getRAbundVector(rank); int numInd = rdata.getNumSeqs(); @@ -65,14 +64,15 @@ EstOutput Geom::getValues(SAbundVector* rank){ double sumExp = 0; double sumObs = 0; double maxDiff = 0; - - for(int i = 0; i < rdata.size(); i++) + + for(int i = 0; i < numSpec; i++) { sumObs += rdata.get(i); sumExp += numInd*cK*k*pow(1-k, i); + double diff = fabs(sumObs-sumExp); - if(diff > maxDiff) - maxDiff = diff; + if(diff > maxDiff) { maxDiff = diff; } + } double DStatistic = maxDiff/numInd; @@ -82,7 +82,8 @@ EstOutput Geom::getValues(SAbundVector* rank){ cout << "Critical value for 95% confidence interval = ";*/ if(rdata.size() > 20) { - critVal = .886/sqrt(rdata.size()); + data[1] = 1.031/sqrt(rdata.size()); + data[2] = 0.886/sqrt(rdata.size()); } else { @@ -93,10 +94,10 @@ EstOutput Geom::getValues(SAbundVector* rank){ cout << "If D-Statistic is less than the critical value then the data fits the Geometric Series model w/ 95% confidence.\n\n";*/ data[0] = DStatistic; - data[1] = critVal; if (isnan(data[0]) || isinf(data[0])) { data[0] = 0; } if (isnan(data[1]) || isinf(data[1])) { data[1] = 0; } + if (isnan(data[2]) || isinf(data[2])) { data[2] = 0; } return data; }