+int CorrAxesCommand::calcSpearman(map<string, vector<float> >& axes, ofstream& out) {
+ try {
+
+ //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++) {
+ vector<float> temp = it->second;
+ for (int i = 0; i < temp.size(); i++) {
+ averageAxes[i] += temp[i];
+ }
+ }
+
+ for (int i = 0; i < averageAxes.size(); i++) { averageAxes[i] = averageAxes[i] / (float) axes.size(); }
+
+ //for each otu
+ for (int i = 0; i < lookupFloat[0]->getNumBins(); i++) {
+
+ out << i+1 << '\t';
+
+ //find the averages this otu - Y
+ float sumOtu = 0.0;
+ for (int j = 0; j < lookupFloat.size(); j++) {
+ sumOtu += lookupFloat[j]->getAbundance(i);
+ }
+ float Ybar = sumOtu / (float) lookupFloat.size();
+
+ //find r value for each axis
+ for (int k = 0; k < averageAxes.size(); k++) {
+
+ double r = 0.0;
+ double numerator = 0.0;
+ double denomTerm1 = 0.0;
+ double denomTerm2 = 0.0;
+ for (int j = 0; j < lookupFloat.size(); j++) {
+ float Yi = lookupFloat[j]->getAbundance(i);
+ float Xi = axes[lookupFloat[j]->getGroup()][k];
+
+ numerator += ((Xi - averageAxes[k]) * (Yi - Ybar));
+ denomTerm1 += ((Xi - averageAxes[k]) * (Xi - averageAxes[k]));
+ denomTerm2 += ((Yi - Ybar) * (Yi - Ybar));
+ }
+
+ double denom = (sqrt(denomTerm1 * denomTerm2));
+
+ r = numerator / denom;
+
+ out << r << '\t';
+ }
+
+ out << endl;
+ }
+
+ return 0;
+ }
+ catch(exception& e) {
+ m->errorOut(e, "CorrAxesCommand", "calcSpearman");
+ exit(1);
+ }
+}
+//**********************************************************************************************************************