X-Git-Url: https://git.donarmstrong.com/?a=blobdiff_plain;f=corraxescommand.cpp;h=c27eb4b9ac335112f19244c708cee9bb8630feda;hb=e7ae6e6b27c45b5691c19f423ec56faae8e2f255;hp=4e9c3c6270b9142ba69963ab0ee2283668473a8b;hpb=0ca63a8165baa0afa459e644ebe140ba496d5ba0;p=mothur.git

diff --git a/corraxescommand.cpp b/corraxescommand.cpp
index 4e9c3c6..c27eb4b 100644
--- a/corraxescommand.cpp
+++ b/corraxescommand.cpp
@@ -9,6 +9,7 @@
 
 #include "corraxescommand.h"
 #include "sharedutilities.h"
+#include "linearalgebra.h"
 
 //**********************************************************************************************************************
 vector<string> CorrAxesCommand::setParameters(){	
@@ -275,7 +276,7 @@ int CorrAxesCommand::execute(){
 		if (metadatafile == "") {  out << "OTU";	}
 		else {  out << "Feature";						}
 
-		for (int i = 0; i < numaxes; i++) { out << '\t' << "axis" << (i+1); }
+		for (int i = 0; i < numaxes; i++) { out << '\t' << "axis" << (i+1) << "\tp-value"; }
 		out << "\tlength" << endl;
 		
 		if (method == "pearson")		{  calcPearson(axes, out);	}
@@ -304,6 +305,8 @@ int CorrAxesCommand::execute(){
 int CorrAxesCommand::calcPearson(map<string, vector<float> >& axes, ofstream& out) {
    try {
 	   
+       LinearAlgebra linear;
+       
 	   //find average of each axis - X
 	   vector<float> averageAxes; averageAxes.resize(numaxes, 0.0);
 	   for (map<string, vector<float> >::iterator it = axes.begin(); it != axes.end(); it++) {
@@ -318,7 +321,7 @@ int CorrAxesCommand::calcPearson(map<string, vector<float> >& axes, ofstream& ou
 	   //for each otu
 	   for (int i = 0; i < lookupFloat[0]->getNumBins(); i++) {
 		   
-		   if (metadatafile == "") {  out << i+1;	}
+		   if (metadatafile == "") {  out << m->currentBinLabels[i];	}
 		   else {  out << metadataLabels[i];		}
 		   		   
 		   //find the averages this otu - Y
@@ -349,8 +352,15 @@ int CorrAxesCommand::calcPearson(map<string, vector<float> >& axes, ofstream& ou
 			   double denom = (sqrt(denomTerm1) * sqrt(denomTerm2));
 			   
 			   r = numerator / denom;
+               
+               if (isnan(r) || isinf(r)) { r = 0.0; }
+               
 			   rValues[k] = r;
 			   out << '\t' << r; 
+               
+               double sig = linear.calcPearsonSig(lookupFloat.size(), r);
+               
+               out << '\t' << sig;
 		   }
 		   
 		   double sum = 0;
@@ -371,6 +381,9 @@ int CorrAxesCommand::calcPearson(map<string, vector<float> >& axes, ofstream& ou
 int CorrAxesCommand::calcSpearman(map<string, vector<float> >& axes, ofstream& out) {
 	try {
 		
+        LinearAlgebra linear;
+        vector<double> sf; 
+        
 		//format data
 		vector< map<float, int> > tableX; tableX.resize(numaxes);
 		map<float, int>::iterator itTable;
@@ -410,6 +423,7 @@ int CorrAxesCommand::calcSpearman(map<string, vector<float> >& axes, ofstream& o
 			
 			vector<spearmanRank> ties;
 			int rankTotal = 0;
+            double sfTemp = 0.0;
 			for (int j = 0; j < scores[i].size(); j++) {
 				rankTotal += (j+1);
 				ties.push_back(scores[i][j]);
@@ -421,6 +435,8 @@ int CorrAxesCommand::calcSpearman(map<string, vector<float> >& axes, ofstream& o
 							float thisrank = rankTotal / (float) ties.size();
   							rankAxes[ties[k].name].push_back(thisrank);
 						}
+                        int t = ties.size();
+                        sfTemp += (t*t*t-t);
 						ties.clear();
 						rankTotal = 0;
 					}
@@ -433,13 +449,14 @@ int CorrAxesCommand::calcSpearman(map<string, vector<float> >& axes, ofstream& o
 					}
 				}
 			}
+            sf.push_back(sfTemp);
 		}
 		
 				
 		//for each otu
 		for (int i = 0; i < lookupFloat[0]->getNumBins(); i++) {
 			
-			if (metadatafile == "") {  out << i+1;	}
+			if (metadatafile == "") {  out << m->currentBinLabels[i];	}
 			else {  out << metadataLabels[i];		}
 			
 			//find the ranks of this otu - Y
@@ -467,6 +484,7 @@ int CorrAxesCommand::calcSpearman(map<string, vector<float> >& axes, ofstream& o
 			
 			sort(otuScores.begin(), otuScores.end(), compareSpearman);
 			
+            double sg = 0.0;
 			map<string, float> rankOtus;
 			vector<spearmanRank> ties;
 			int rankTotal = 0;
@@ -481,6 +499,8 @@ int CorrAxesCommand::calcSpearman(map<string, vector<float> >& axes, ofstream& o
 							float thisrank = rankTotal / (float) ties.size();
   							rankOtus[ties[k].name] = thisrank;
 						}
+                        int t = ties.size();
+                        sg += (t*t*t-t);
 						ties.clear();
 						rankTotal = 0;
 					}
@@ -515,9 +535,15 @@ int CorrAxesCommand::calcSpearman(map<string, vector<float> >& axes, ofstream& o
 				
 				p = (SX2 + SY2 - di) / (2.0 * sqrt((SX2*SY2)));
 				
+                if (isnan(p) || isinf(p)) { p = 0.0; }
+                
 				out  << '\t' << p;
 				
 				pValues[j] = p;
+                
+                double sig = linear.calcSpearmanSig(n, sf[j], sg, di);            
+                out  << '\t' << sig;
+                
 			}
 
 			double sum = 0;
@@ -538,6 +564,8 @@ int CorrAxesCommand::calcSpearman(map<string, vector<float> >& axes, ofstream& o
 int CorrAxesCommand::calcKendall(map<string, vector<float> >& axes, ofstream& out) {
 	try {
 		
+        LinearAlgebra linear;
+        
 		//format data
 		vector< vector<spearmanRank> > scores; scores.resize(numaxes);
 		for (map<string, vector<float> >::iterator it = axes.begin(); it != axes.end(); it++) {
@@ -581,7 +609,7 @@ int CorrAxesCommand::calcKendall(map<string, vector<float> >& axes, ofstream& ou
 		//for each otu
 		for (int i = 0; i < lookupFloat[0]->getNumBins(); i++) {
 		
-			if (metadatafile == "") {  out << i+1;	}
+			if (metadatafile == "") {  out << m->currentBinLabels[i];	}
 			else {  out << metadataLabels[i];		}
 			
 			//find the ranks of this otu - Y
@@ -651,10 +679,14 @@ int CorrAxesCommand::calcKendall(map<string, vector<float> >& axes, ofstream& ou
 				}
 				
 				double p = (numCoor - numDisCoor) / (float) count;
-
+                 if (isnan(p) || isinf(p)) { p = 0.0; }
+                
 				out << '\t' << p;
 				pValues[j] = p;
-
+                
+                double sig = linear.calcKendallSig(scores[j].size(), p);
+                
+                out << '\t' << sig;
 			}
 			
 			double sum = 0;