changes for 1.16.0

diff --git a/Mothur.xcodeproj/project.pbxproj b/Mothur.xcodeproj/project.pbxproj
index 112cf3c7fdc7305c62bce17667f81232492860f6..c1a1dc9e612f705594a1538247f4586873a7e931 100644 (file)
--- a/Mothur.xcodeproj/project.pbxproj
+++ b/Mothur.xcodeproj/project.pbxproj
@@ -958,6 +958,12 @@
                                A7E9B75C12D37EC400DA6239 /* mothur.h */,
                                A7E9B75D12D37EC400DA6239 /* mothurout.cpp */,
                                A7E9B75E12D37EC400DA6239 /* mothurout.h */,
+                               A7E9B76112D37EC400DA6239 /* nast.cpp */,
+                               A7E9B76212D37EC400DA6239 /* nast.hpp */,
+                               A7E9B76312D37EC400DA6239 /* nastreport.cpp */,
+                               A7E9B76412D37EC400DA6239 /* nastreport.hpp */,
+                               A7E9B76712D37EC400DA6239 /* noalign.cpp */,
+                               A7E9B76812D37EC400DA6239 /* noalign.hpp */,
                                A7E9B76512D37EC400DA6239 /* needlemanoverlap.cpp */,
                                A7E9B76612D37EC400DA6239 /* needlemanoverlap.hpp */,
                                A7E9B77012D37EC400DA6239 /* observable.h */,
@@ -974,13 +980,9 @@
                                A7E9B7A912D37EC400DA6239 /* rarefact.cpp */,
                                A7E9B7AA12D37EC400DA6239 /* rarefact.h */,
                                A7E9B7AD12D37EC400DA6239 /* rarefactioncurvedata.h */,
+                               7E6BE10812F710D8007ADDBE /* refchimeratest.h */,
+                               7E6BE10912F710D8007ADDBE /* refchimeratest.cpp */,
                                A7E9BA5312D39A5E00DA6239 /* read */,
-                               A7E9B76112D37EC400DA6239 /* nast.cpp */,
-                               A7E9B76212D37EC400DA6239 /* nast.hpp */,
-                               A7E9B76312D37EC400DA6239 /* nastreport.cpp */,
-                               A7E9B76412D37EC400DA6239 /* nastreport.hpp */,
-                               A7E9B76712D37EC400DA6239 /* noalign.cpp */,
-                               A7E9B76812D37EC400DA6239 /* noalign.hpp */,
                                A7E9B82312D37EC400DA6239 /* sharedutilities.cpp */,
                                A7E9B82412D37EC400DA6239 /* sharedutilities.h */,
                                A7E9B82D12D37EC400DA6239 /* singlelinkage.cpp */,
@@ -1014,8 +1016,6 @@
                A7E9BA3812D3956100DA6239 /* commands */ = {
                        isa = PBXGroup;
                        children = (
-                               7E6BE10812F710D8007ADDBE /* refchimeratest.h */,
-                               7E6BE10912F710D8007ADDBE /* refchimeratest.cpp */,
                                A7E9B6AE12D37EC400DA6239 /* command.hpp */,
                                A7E9B65112D37EC300DA6239 /* aligncommand.cpp */,
                                A7E9B65212D37EC300DA6239 /* aligncommand.h */,
@@ -1598,7 +1598,7 @@
                        attributes = {
                                ORGANIZATIONNAME = "Schloss Lab";
                        };
-                       buildConfigurationList = 1DEB928908733DD80010E9CD /* Build configuration list for PBXProject "Mothur" */;
+                       buildConfigurationList = 1DEB928908733DD80010E9CD /* Build configuration list for PBXProject "mothur" */;
                        compatibilityVersion = "Xcode 3.1";
                        developmentRegion = English;
                        hasScannedForEncodings = 1;
@@ -2016,7 +2016,7 @@
                        defaultConfigurationIsVisible = 0;
                        defaultConfigurationName = Release;
                };
-               1DEB928908733DD80010E9CD /* Build configuration list for PBXProject "Mothur" */ = {
+               1DEB928908733DD80010E9CD /* Build configuration list for PBXProject "mothur" */ = {
                        isa = XCConfigurationList;
                        buildConfigurations = (
                                1DEB928A08733DD80010E9CD /* Debug */,

diff --git a/catchallcommand.cpp b/catchallcommand.cpp
index 8e29d6294c942369b178eea7b3cf6a0c3e7ba272..69b9f04becdbb492fc69209837be5a0581123b99 100644 (file)
--- a/catchallcommand.cpp
+++ b/catchallcommand.cpp
@@ -121,6 +121,11 @@ CatchAllCommand::CatchAllCommand(string option)  {
                        if (sharedfile == "not open") { sharedfile = ""; abort = true; }
                        else if (sharedfile == "not found") { sharedfile = "";   }
                        
+                       //check for shared file loaded during read.otu
+                       if (sharedfile == "") {
+                               if (globaldata->getSharedFile() != "") { sharedfile = globaldata->getSharedFile(); }
+                       }
+                       
                        string label = validParameter.validFile(parameters, "label", false);                    
                        if (label == "not found") { label = ""; }
                        else { 

diff --git a/corraxescommand.cpp b/corraxescommand.cpp
index db2e9891738cdf57b78dffe1a9c232f18e5de31c..eac1f3505b1bca098ab01d100a513d24bd5c9d29 100644 (file)
--- a/corraxescommand.cpp
+++ b/corraxescommand.cpp
@@ -382,6 +382,8 @@ int CorrAxesCommand::calcPearson(map<string, vector<float> >& axes, ofstream& ou
 int CorrAxesCommand::calcSpearman(map<string, vector<float> >& axes, ofstream& out) {
        try {
                
+               bool hasTies = false;
+               
                //format data
                vector< vector<spearmanRank> > scores; scores.resize(numaxes);
                for (map<string, vector<float> >::iterator it = axes.begin(); it != axes.end(); it++) {
@@ -396,6 +398,7 @@ int CorrAxesCommand::calcSpearman(map<string, vector<float> >& axes, ofstream& o
                for (int i = 0; i < numaxes; i++) {  sort(scores[i].begin(), scores[i].end(), compareSpearman); }
                
                //find ranks of xi in each axis
+               vector<float> averageX; averageX.resize(numaxes, 0.0);
                map<string, vector<float> > rankAxes;
                for (int i = 0; i < numaxes; i++) {
                        
@@ -407,6 +410,8 @@ int CorrAxesCommand::calcSpearman(map<string, vector<float> >& axes, ofstream& o
                                
                                if (j != scores.size()-1) { // you are not the last so you can look ahead
                                        if (scores[i][j].score != scores[i][j+1].score) { // you are done with ties, rank them and continue
+                                               if (ties.size() > 1) { hasTies = true; }
+                                               averageX[i] += rankTotal;
                                                for (int k = 0; k < ties.size(); k++) {
                                                        float thisrank = rankTotal / (float) ties.size();
                                                        rankAxes[ties[k].name].push_back(thisrank);
@@ -415,12 +420,16 @@ int CorrAxesCommand::calcSpearman(map<string, vector<float> >& axes, ofstream& o
                                                rankTotal = 0;
                                        }
                                }else { // you are the last one
+                                       if (ties.size() > 1) { hasTies = true; }
+                                       averageX[i] += rankTotal;
                                        for (int k = 0; k < ties.size(); k++) {
                                                float thisrank = rankTotal / (float) ties.size();
                                                rankAxes[ties[k].name].push_back(thisrank);
                                        }
                                }
                        }
+                       
+                       averageX[i] /= (float) scores[i].size();
                }
                
                                
@@ -442,6 +451,7 @@ int CorrAxesCommand::calcSpearman(map<string, vector<float> >& axes, ofstream& o
                        
                        map<string, float> rankOtus;
                        vector<spearmanRank> ties;
+                       float averageY = 0.0;
                        int rankTotal = 0;
                        for (int j = 0; j < otuScores.size(); j++) {
                                rankTotal += j;
@@ -449,6 +459,8 @@ int CorrAxesCommand::calcSpearman(map<string, vector<float> >& axes, ofstream& o
                                
                                if (j != scores.size()-1) { // you are not the last so you can look ahead
                                        if (otuScores[j].score != otuScores[j+1].score) { // you are done with ties, rank them and continue
+                                               if (ties.size() > 1) { hasTies = true; }
+                                               averageY += rankTotal;
                                                for (int k = 0; k < ties.size(); k++) {
                                                        float thisrank = rankTotal / (float) ties.size();
                                                        rankOtus[ties[k].name] = thisrank;
@@ -457,6 +469,8 @@ int CorrAxesCommand::calcSpearman(map<string, vector<float> >& axes, ofstream& o
                                                rankTotal = 0;
                                        }
                                }else { // you are the last one
+                                       if (ties.size() > 1) { hasTies = true; }
+                                       averageY += rankTotal;
                                        for (int k = 0; k < ties.size(); k++) {
                                                float thisrank = rankTotal / (float) ties.size();
                                                rankOtus[ties[k].name] = thisrank;
@@ -464,32 +478,44 @@ int CorrAxesCommand::calcSpearman(map<string, vector<float> >& axes, ofstream& o
                                }
                        }
                        
-                       vector<double> pValues(numaxes);
+
+                       averageY /= (float) otuScores.size();
+                       //hasTies = false;              
+
                        //calc spearman ranks for each axis for this otu
                        for (int j = 0; j < numaxes; j++) {
                                
                                double di = 0.0;
+                               double denom1 = 0.0;
+                               double denom2 = 0.0;
                                for (int k = 0; k < lookupFloat.size(); k++) {
                                        
                                        float xi = rankAxes[lookupFloat[k]->getGroup()][j];
                                        float yi = rankOtus[lookupFloat[k]->getGroup()];
                                        
-                                       di += ((xi - yi) * (xi - yi));
+                                       if (hasTies) {
+                                               di += ((xi - averageX[j]) * (yi - averageY));
+                                               denom1 += ((xi - averageX[j]) * (xi - averageX[j]));
+                                               denom2 += ((yi - averageY) * (yi - averageY));
+                                       }else {
+                                               di += ((xi - yi) * (xi - yi));
+                                       }
                                }
                                
-                               int n = lookupFloat.size();
-                               double p = 1.0 - ((6 * di) / (float) (n * ((n*n) - 1)));
+                               double p = 0.0;
+                               
+                               if (hasTies) {
+                                       double denom = sqrt((denom1 * denom2));
+                                       p = di / denom;
+                               }else {
+                                       int n = lookupFloat.size();
+                                       p = 1.0 - ((6 * di) / (float) (n * ((n*n) - 1)));
+                               }
                                
                                out  << '\t' << p;
-                               pValues[j] = p;
-
                        }
 
-                       double sum = 0;
-                       for(int k=0;k<numaxes;k++){
-                               sum += pValues[k] * pValues[k];
-                       }
-                       out << '\t' << sqrt(sum) << endl;
+                       out << endl;
                }
                
                return 0;

diff --git a/splitabundcommand.cpp b/splitabundcommand.cpp
index fb0326a479efbc64d09bf946316355ce8dc44bac..80145bfee4f3bfad04ae5becb2f5198e7a109abb 100644 (file)
--- a/splitabundcommand.cpp
+++ b/splitabundcommand.cpp
@@ -237,28 +237,6 @@ int SplitAbundCommand::execute(){
                
                if (listfile != "") { //you are using a listfile to determine abundance
                        if (outputDir == "") { outputDir = m->hasPath(listfile); }
-               
-                       //remove old files so you can append later....
-                       string fileroot = outputDir + m->getRootName(m->getSimpleName(listfile));
-                       if (Groups.size() == 0) {
-                               remove((fileroot + "rare.list").c_str());
-                               remove((fileroot + "abund.list").c_str());
-                               
-                               outputNames.push_back((fileroot + "rare.list"));
-                               outputNames.push_back((fileroot + "abund.list"));
-                               outputTypes["list"].push_back((fileroot + "rare.list"));
-                               outputTypes["list"].push_back((fileroot + "abund.list"));
-                       }else{
-                               for (int i=0; i<Groups.size(); i++) {
-                                       remove((fileroot + Groups[i] + ".rare.list").c_str());
-                                       remove((fileroot + Groups[i] + ".abund.list").c_str());
-                                       
-                                       outputNames.push_back((fileroot + Groups[i] + ".rare.list"));
-                                       outputNames.push_back((fileroot + Groups[i] + ".abund.list"));
-                                       outputTypes["list"].push_back((fileroot + Groups[i] + ".rare.list"));
-                                       outputTypes["list"].push_back((fileroot + Groups[i] + ".abund.list"));
-                               }
-                       }
                        
                        //if the users enters label "0.06" and there is no "0.06" in their file use the next lowest label.
                        set<string> processedLabels;
@@ -390,9 +368,11 @@ int SplitAbundCommand::splitList(ListVector* thisList) {
                        }
                }//end for
 
-               writeList(thisList);
                
                string tag = thisList->getLabel() + ".";
+               
+               writeList(thisList, tag);
+               
                if (groupfile != "")                            {  parseGroup(tag);             }
                if (accnos)                                                     {  writeAccnos(tag);    }
                if (fastafile != "")                            {  parseFasta(tag);             }
@@ -406,7 +386,7 @@ int SplitAbundCommand::splitList(ListVector* thisList) {
        }
 }
 /**********************************************************************************************************************/
-int SplitAbundCommand::writeList(ListVector* thisList) { 
+int SplitAbundCommand::writeList(ListVector* thisList, string tag) { 
        try {
                
                map<string, ofstream*> filehandles;
@@ -428,13 +408,13 @@ int SplitAbundCommand::writeList(ListVector* thisList) {
                        ofstream aout;
                        ofstream rout;
                        
-                       string rare = outputDir + m->getRootName(m->getSimpleName(listfile))  + "rare.list";
-                       m->openOutputFileAppend(rare, rout);
-                       //outputNames.push_back(rare);
+                       string rare = outputDir + m->getRootName(m->getSimpleName(listfile)) + tag + "rare.list";
+                       m->openOutputFile(rare, rout);
+                       outputNames.push_back(rare); outputTypes["list"].push_back(rare);
                        
-                       string abund = outputDir + m->getRootName(m->getSimpleName(listfile))  + "abund.list";
-                       m->openOutputFileAppend(abund, aout);
-                       //outputNames.push_back(abund);
+                       string abund = outputDir + m->getRootName(m->getSimpleName(listfile)) + tag + "abund.list";
+                       m->openOutputFile(abund, aout);
+                       outputNames.push_back(abund); outputTypes["list"].push_back(abund);
 
                        if (rareNames.size() != 0)      {  rout << thisList->getLabel() << '\t' << numRareBins << '\t';         }
                        if (abundNames.size() != 0) {   aout << thisList->getLabel() << '\t' << numAbundBins << '\t';   }
@@ -470,8 +450,10 @@ int SplitAbundCommand::writeList(ListVector* thisList) {
                                temp2 = new ofstream;
                                filehandles[Groups[i]+".abund"] = temp2;
                                
-                               m->openOutputFileAppend(fileroot + Groups[i] + ".rare.list", *(filehandles[Groups[i]+".rare"]));
-                               m->openOutputFileAppend(fileroot + Groups[i] + ".abund.list", *(filehandles[Groups[i]+".abund"]));
+                               m->openOutputFile(fileroot + Groups[i] + tag + ".rare.list", *(filehandles[Groups[i]+".rare"]));
+                               m->openOutputFile(fileroot + Groups[i] + tag + ".abund.list", *(filehandles[Groups[i]+".abund"]));
+                               outputNames.push_back(fileroot + Groups[i] + tag + ".rare.list"); outputTypes["list"].push_back(fileroot + Groups[i] + tag + ".rare.list");
+                               outputNames.push_back(fileroot + Groups[i] + tag + ".abund.list"); outputTypes["list"].push_back(fileroot + Groups[i] + tag + ".abund.list");
                        }
                        
                        map<string, string> groupVector;

diff --git a/splitabundcommand.h b/splitabundcommand.h
index 02d61f3dc3b7f343a9e0d9a73bb71e84fcdacbbd..2e409b2aedd7606b6c689f48cb9e0e58f7a0bcfe 100644 (file)
--- a/splitabundcommand.h
+++ b/splitabundcommand.h
@@ -42,7 +42,7 @@ private:
        int splitList(ListVector*);
        int splitNames(); //namefile
        int writeNames(); 
-       int writeList(ListVector*); 
+       int writeList(ListVector*, string); 
        int writeAccnos(string); 
        int parseGroup(string); 
        int parseFasta(string);