5 * Created by Sarah Westcott on 3/9/09.
6 * Copyright 2009 Schloss Lab UMASS Amherst. All rights reserved.
10 /* This class library coverage at the given distances of the Cramer-von Mises statistic.
11 you may refer to the "Integration of Microbial Ecology and Statistics: A Test To Compare Gene Libraries"
12 paper in Applied and Environmental Microbiology, Sept. 2004, p. 5485-5492 0099-2240/04/$8.00+0
13 DOI: 10.1128/AEM.70.9.5485-5492.2004 Copyright 2004 American Society for Microbiology for more information. */
18 //**********************************************************************************************************************
19 Coverage::Coverage() {
20 globaldata = GlobalData::getInstance();
21 numUserGroups = globaldata->Groups.size();
22 numGroups = globaldata->gGroupmap->getNumGroups();
25 //**********************************************************************************************************************
26 void Coverage::getValues(FullMatrix* matrix, vector< vector< vector<float> > >& data, vector<float> dist) {
29 vector<string> groups;
32 data.resize(dist.size());
33 for (int l = 0; l < data.size(); l++) {
34 data[l].resize(numGroups);
35 for (int k = 0; k < data[l].size(); k++) {
36 data[l][k].push_back(0.0);
40 /**************************************/
41 //get the minimums for each comparision
42 /**************************************/
45 for (int i = 0; i < numGroups; i++) {
46 for (int j = 0; j < numGroups; j++) {
48 //is this "box" one the user wants analyzed?
49 if ((inUsersGroups(globaldata->gGroupmap->namesOfGroups[i], globaldata->Groups) == true) && (inUsersGroups(globaldata->gGroupmap->namesOfGroups[j], globaldata->Groups) == true)) {
51 min = matrix->getMins(count); //returns vector of mins for "box" requested ie. groups A, B, 0 = AA, 1 = AB, 2 = BA, 3 = BB;
53 //find the coverage at this distance
54 sort(min.begin(), min.end());
56 //loop through each distance and fill data
57 for (int k = 0; k < data.size(); k++) {
60 //find index in min where value is higher than d
61 for (int m = 0; m < min.size(); m++) {
62 if (min[m] > dist[k]) { index = m; break; }
65 // if you don't find one than all the mins are less than d
66 if (index == -1) { index = min.size(); }
69 data[k][i][j] = 1.0 - ((min.size()-index)/(float)min.size());
79 cout << "Standard Error: " << e.what() << " has occurred in the Coverage class Function getValues. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
83 cout << "An unknown error has occurred in the Coverage class function getValues. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
88 //**********************************************************************************************************************
89 void Coverage::getValues(FullMatrix* matrix, vector< vector< vector<float> > >& data, vector<float> dist, string mode) {
93 vector<string> groups;
96 data.resize(dist.size());
97 for (int l = 0; l < data.size(); l++) {
98 data[l].resize(numGroups);
99 for (int k = 0; k < data[l].size(); k++) {
100 data[l][k].push_back(0.0);
104 /**************************************/
105 //get the minimums for each comparision
106 /**************************************/
111 for (int i = 0; i < numGroups; i++) {
112 for (int j = 0; j < numGroups; j++) {
115 //is this "box" one the user wants analyzed?
116 if ((inUsersGroups(globaldata->gGroupmap->namesOfGroups[i], globaldata->Groups) == true) && (inUsersGroups(globaldata->gGroupmap->namesOfGroups[j], globaldata->Groups) == true)) {
118 matrix->shuffle(globaldata->gGroupmap->namesOfGroups[i], globaldata->gGroupmap->namesOfGroups[j]);
120 min1 = matrix->getMins(count); //returns vector of mins for "box" requested ie. groups A, B, 0 = AA, 1 = AB, 2 = BA, 3 = BB;
121 min2 = matrix->getMins(count2); //returns vector of mins for "box" requested ie. groups A, B, 0 = AA, 1 = AB, 2 = BA, 3 = BB;
123 //find the coverage at this distance
124 sort(min1.begin(), min1.end());
126 //find the coverage at this distance
127 sort(min2.begin(), min2.end());
131 //loop through each distance and fill data
132 for (int k = 0; k < data.size(); k++) {
133 //****** coverage of AA **********//
135 //find index in min where value is higher than d
136 for (int m = 0; m < min1.size(); m++) {
137 if (min1[m] > dist[k]) { index = m; break; }
140 // if you don't find one than all the mins are less than d
141 if (index == -1) { index = min1.size(); }
143 //****** coverage of AB **********//
145 //find index in min where value is higher than d
146 for (int m = 0; m < min2.size(); m++) {
147 if (min2[m] > dist[k]) { index2 = m; break; }
150 // if you don't find one than all the mins are less than d
151 if (index2 == -1) { index2 = min2.size(); }
154 float covII = 1.0 - ((min1.size()-index)/(float)min1.size());
157 float covIJ = 1.0 - ((min2.size()-index2)/(float)min2.size());
159 //save value in data (Caa - Cab)^2 * distDiff
160 data[k][i][j] = ((covII-covIJ) * (covII-covIJ)) * distDiff;
163 if (k < data.size() - 1) {
164 distDiff = dist[k+1] - dist[k];
168 //put matrix back to original
175 count += numGroups+1; //go from AA to BB to CC
179 catch(exception& e) {
180 cout << "Standard Error: " << e.what() << " has occurred in the Coverage class Function getValues. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
184 cout << "An unknown error has occurred in the Coverage class function getValues. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";