From: westcott <westcott>
Date: Thu, 3 Feb 2011 16:44:19 +0000 (+0000)
Subject: changes for 1.16.0
X-Git-Url: https://git.donarmstrong.com/?a=commitdiff_plain;h=0e40e23448c2ee274268d85e0d0e65cb9eaeee6f;p=mothur.git

changes for 1.16.0
---

diff --git a/Mothur.xcodeproj/project.pbxproj b/Mothur.xcodeproj/project.pbxproj
index 112cf3c..c1a1dc9 100644
--- 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 8e29d62..69b9f04 100644
--- 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 db2e989..eac1f35 100644
--- 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 fb0326a..80145bf 100644
--- 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 02d61f3..2e409b2 100644
--- 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);