]> git.donarmstrong.com Git - mothur.git/blob - sharedutilities.cpp
sped up unifrac unweighted.
[mothur.git] / sharedutilities.cpp
1 /*
2  *  sharedutilities.cpp
3  *  Mothur
4  *
5  *  Created by Sarah Westcott on 4/9/09.
6  *  Copyright 2009 Schloss Lab UMASS Amherst. All rights reserved.
7  *
8  */
9
10 #include "sharedutilities.h"
11 #include "sharedrabundvector.h"
12 #include "sharedordervector.h"
13
14 /**************************************************************************************************/
15
16 void SharedUtil::getSharedVectors(vector<string> Groups, vector<SharedRAbundVector*>& lookup, SharedOrderVector* order) {
17         try {
18         
19                 //delete each sharedrabundvector in lookup
20                 for (int j = 0; j < lookup.size(); j++) {
21                         delete lookup[j];
22                 }
23                 
24                 lookup.clear();
25                 
26                 //create and initialize vector of sharedvectors, one for each group
27                 for (int i = 0; i < Groups.size(); i++) { 
28                         SharedRAbundVector* temp = new SharedRAbundVector(order->getNumBins());
29                         temp->setLabel(order->getLabel());
30                         temp->setGroup(Groups[i]);
31                         lookup.push_back(temp);
32                 }
33         
34                 int numSeqs = order->size();
35                 //sample all the members
36                 for(int i=0;i<numSeqs;i++){
37                         //get first sample
38                         individual chosen = order->get(i);
39                         int abundance; 
40                                         
41                         //set info for sharedvector in chosens group
42                         for (int j = 0; j < lookup.size(); j++) { 
43                                 if (chosen.group == lookup[j]->getGroup()) {
44                                          abundance = lookup[j]->getAbundance(chosen.bin);
45                                          lookup[j]->set(chosen.bin, (abundance + 1), chosen.group);
46                                          break;
47                                 }
48                         }
49                 }
50         }
51         catch(exception& e) {
52                 m->errorOut(e, "SharedUtil", "getSharedVectors");
53                 exit(1);
54         }
55 }
56 /**************************************************************************************************/
57
58 void SharedUtil::getSharedVectorswithReplacement(vector<string> Groups, vector<SharedRAbundVector*>& lookup, SharedOrderVector* order) {
59         try {
60         
61                 //delete each sharedrabundvector in lookup
62                 for (int j = 0; j < lookup.size(); j++) {
63                         delete lookup[j];
64                 }
65                 lookup.clear();
66                 
67                 //create and initialize vector of sharedvectors, one for each group
68                 for (int i = 0; i < Groups.size(); i++) { 
69                         SharedRAbundVector* temp = new SharedRAbundVector(order->getNumBins());
70                         temp->setLabel(order->getLabel());
71                         temp->setGroup(Groups[i]);
72                         lookup.push_back(temp);
73                 }
74         
75                 int numSeqs = order->size();
76                 
77                 //sample all the members
78                 for(int i=0;i<numSeqs;i++){
79                         //get random number
80                         int random = int((float)(i+1) * (float)(rand()) / ((float)RAND_MAX+1.0));
81                         individual chosen = order->get(random);
82
83                         int abundance; 
84                         //set info for sharedvector in chosens group
85                         for (int j = 0; j < lookup.size(); j++) { 
86                                 if (chosen.group == lookup[j]->getGroup()) {
87                                          abundance = lookup[j]->getAbundance(chosen.bin);
88                                          lookup[j]->set(chosen.bin, (abundance + 1), chosen.group);
89                                          break;
90                                 }
91                         }
92                 }
93                 
94         }
95         catch(exception& e) {
96                 m->errorOut(e, "SharedUtil", "getSharedVectorswithReplacement");
97                 exit(1);
98         }
99 }
100
101 /**************************************************************************************************/
102 //need to have mode because different commands require different number of valid groups
103 void SharedUtil::setGroups(vector<string>& userGroups, vector<string>& allGroups) {
104         try {
105                 if (userGroups.size() != 0) {
106                         if (userGroups[0] != "all") {
107                                 //check that groups are valid
108                                 for (int i = 0; i < userGroups.size(); i++) {
109                                         if (isValidGroup(userGroups[i], allGroups) != true) {
110                                                 m->mothurOut(userGroups[i] + " is not a valid group, and will be disregarded."); m->mothurOutEndLine();
111                                                 // erase the invalid group from userGroups
112                                                 userGroups.erase(userGroups.begin()+i);
113                                                 i--;
114                                         }
115                                 }
116                                 
117                                 //if the user only entered invalid groups
118                                 if (userGroups.size() == 0) { 
119                                         m->mothurOut("You provided no valid groups. I will run the command using all the groups in your groupfile."); m->mothurOutEndLine();
120                                         for (int i = 0; i < allGroups.size(); i++) {
121                                                 userGroups.push_back(allGroups[i]);
122                                         }
123                                 }
124
125                         }else{//user has enter "all" and wants the default groups
126                                 userGroups.clear();
127                                 for (int i = 0; i < allGroups.size(); i++) {
128                                         userGroups.push_back(allGroups[i]);
129                                 }
130                         }
131                 }else { //the user has not entered groups
132                         for (int i = 0; i < allGroups.size(); i++) {
133                                 userGroups.push_back(allGroups[i]);
134                         }
135                 }
136                         
137         }
138         catch(exception& e) {
139                 m->errorOut(e, "SharedUtil", "setGroups");
140                 exit(1);
141         }
142 }
143 /**************************************************************************************************/
144 //need to have mode because different commands require different number of valid groups
145 void SharedUtil::setGroups(vector<string>& userGroups, vector<string>& allGroups, string mode) {
146         try {
147                 if (userGroups.size() != 0) {
148                         if (userGroups[0] != "all") {
149                                 //check that groups are valid
150                                 for (int i = 0; i < userGroups.size(); i++) {
151                                         if (isValidGroup(userGroups[i], allGroups) != true) {
152                                                 m->mothurOut(userGroups[i] + " is not a valid group, and will be disregarded."); m->mothurOutEndLine();
153                                                 // erase the invalid group from userGroups
154                                                 userGroups.erase(userGroups.begin()+i);
155                                                 i--;
156                                         }
157                                 }
158
159                         }else{//user has enter "all" and wants the default groups
160                                 userGroups.clear();
161                                 for (int i = 0; i < allGroups.size(); i++) {
162                                         userGroups.push_back(allGroups[i]);
163                                 }
164                         }
165                 }else { //the user has not entered groups
166                         for (int i = 0; i < allGroups.size(); i++) {
167                                 userGroups.push_back(allGroups[i]);
168                         }
169                 }
170                         
171                 if ((mode == "collect") || (mode == "rarefact") || (mode == "summary") || (mode == "treegroup")) {
172                                 //if the user only entered invalid groups
173                                 if ((userGroups.size() == 0) || (userGroups.size() == 1)) { 
174                                         m->mothurOut("When using the groups parameter you must have at least 2 valid groups. I will run the command using all the groups in your groupfile."); m->mothurOutEndLine();
175                                         for (int i = 0; i < allGroups.size(); i++) {
176                                                 userGroups.push_back(allGroups[i]);
177                                         }
178                                 }
179                 }
180         
181         }
182         catch(exception& e) {
183                 m->errorOut(e, "SharedUtil", "setGroups");
184                 exit(1);
185         }
186 }
187
188
189 /**************************************************************************************/
190 //for parsimony and unifrac commands you set pairwise groups as well as an allgroups in calc
191 void SharedUtil::setGroups(vector<string>& userGroups, vector<string>& allGroups, string& label, int& numGroups, string mode){  //globaldata->Groups, your tree or group map, allgroups, mode
192         try {
193                 numGroups = 0;
194                 label = "";
195
196                 //if the user has not entered specific groups to analyze then do them all
197                 if (userGroups.size() != 0) {
198                         if (userGroups[0] != "all") {
199                                 //check that groups are valid
200                                 for (int i = 0; i < userGroups.size(); i++) {
201                                         if (isValidGroup(userGroups[i], allGroups) != true) {
202                                                 m->mothurOut(userGroups[i] + " is not a valid group, and will be disregarded."); m->mothurOutEndLine();
203                                                 // erase the invalid group from globaldata->Groups
204                                                 userGroups.erase(userGroups.begin()+i);
205                                                 i--;
206                                         }
207                                 }
208                         }else { //users wants all groups
209                                 userGroups.clear();
210                                 for (int i=0; i < allGroups.size(); i++) { 
211                                         if (allGroups[i] != "xxx") {
212                                                 userGroups.push_back(allGroups[i]);
213                                         }
214                                 }
215                         }
216                 }else { //the user has not entered groups
217                         for (int i=0; i < allGroups.size(); i++) { 
218                                 if (allGroups[i] != "xxx") {
219                                         if (mode == "weighted") {
220                                                 userGroups.push_back(allGroups[i]);
221                                         }else {
222                                                 numGroups = 1;
223                                                 label += allGroups[i] + "-";
224                                         }
225                                 }
226                         }
227                         //rip extra - off allgroups 
228                         label = label.substr(0, label.length()-1);
229                         if ((mode != "weighted") && (allGroups.size() > 10)) {  label = "merged";  }
230                 }
231                 
232                 if (mode == "weighted") {
233                         //if the user only entered invalid groups
234                         if (userGroups.size() == 0) { 
235                                 for (int i=0; i < allGroups.size(); i++) { 
236                                         if (allGroups[i] != "xxx") {
237                                                 userGroups.push_back(allGroups[i]);
238                                         }
239                                 }
240                                 m->mothurOut("When using the groups parameter you must have at least 2 valid groups. I will run the command using all the groups in your groupfile."); m->mothurOutEndLine();
241                         }else if (userGroups.size() == 1) { 
242                                 m->mothurOut("When using the groups parameter you must have at least 2 valid groups. I will run the command using all the groups in your groupfile."); m->mothurOutEndLine();
243                                 userGroups.clear();
244                                 for (int i=0; i < allGroups.size(); i++) { 
245                                         if (allGroups[i] != "xxx") {
246                                                 userGroups.push_back(allGroups[i]);
247                                         }
248                                 }
249                         }
250                         numGroups = userGroups.size();
251                         
252                 }else if ((mode == "unweighted") || (mode == "parsimony")) {
253                                 //if the user only entered invalid groups
254                                 if ((userGroups.size() == 0) && (numGroups == 0)) { 
255                                         m->mothurOut("When using the groups parameter you must have at least 1 valid group. I will run the command using all the groups in your groupfile."); m->mothurOutEndLine();
256                                         for (int i = 0; i < allGroups.size(); i++) {
257                                                 if (allGroups[i] != "xxx") {
258                                                         userGroups.push_back(allGroups[i]);
259                                                 }
260                                         }
261                                 }
262                                 
263                                 if (numGroups != 1) { numGroups = userGroups.size(); }
264                 }
265         }
266         catch(exception& e) {
267                 m->errorOut(e, "SharedUtil", "setGroups");
268                 exit(1);
269         }
270 }
271 /**************************************************************************************/
272 void SharedUtil::getCombos(vector<string>& groupComb, vector<string> userGroups, int& numComp) { //groupcomb, globaldata->Groups, numcomb
273         try {
274                 //calculate number of comparisons i.e. with groups A,B,C = AB, AC, BC = 3;
275                 numComp = 0;
276                 for (int i=0; i< userGroups.size(); i++) { 
277                         numComp += i; 
278                         for (int l = i+1; l < userGroups.size(); l++) {
279                                 //set group comparison labels
280                                 groupComb.push_back(userGroups[i] + "-" + userGroups[l]);
281                         }
282                 } 
283         }
284         catch(exception& e) {
285                 m->errorOut(e, "SharedUtil", "getCombos");
286                 exit(1);
287         }
288 }
289 /**************************************************************************************/
290 bool SharedUtil::isValidGroup(string groupname, vector<string> groups) {
291         try {
292                 for (int i = 0; i < groups.size(); i++) {
293                         if (groupname == groups[i]) { return true; }
294                 }
295                 
296                 return false;
297         }
298         catch(exception& e) {
299                 m->errorOut(e, "SharedUtil", "isValidGroup");
300                 exit(1);
301         }
302 }
303
304 /**************************************************************************************/
305 void SharedUtil::updateGroupIndex(vector<string>& userGroups, map<string, int>& index) {
306         try {
307                 index.clear();
308                 for (int i = 0; i < userGroups.size(); i++) {
309                         index[userGroups[i]] = i;
310                 }
311         }
312         catch(exception& e) {
313                 m->errorOut(e, "SharedUtil", "updateGroupIndex");
314                 exit(1);
315         }
316 }