//**********************************************************************************************************************
vector<string> IndicatorCommand::setParameters(){
try {
- CommandParameter piters("iters", "Number", "", "1000", "", "", "",false,false); parameters.push_back(piters);
- CommandParameter pdesign("design", "InputTypes", "", "", "TreeDesign", "TreeDesign", "none",false,false); parameters.push_back(pdesign);
- CommandParameter pshared("shared", "InputTypes", "", "", "SharedRel", "SharedRel", "none",false,false); parameters.push_back(pshared);
- CommandParameter prelabund("relabund", "InputTypes", "", "", "SharedRel", "SharedRel", "none",false,false); parameters.push_back(prelabund);
- CommandParameter pgroups("groups", "String", "", "", "", "", "",false,false); parameters.push_back(pgroups);
- CommandParameter plabel("label", "String", "", "", "", "", "",false,false); parameters.push_back(plabel);
- CommandParameter ptree("tree", "InputTypes", "", "", "TreeDesign", "TreeDesign", "none",false,false); parameters.push_back(ptree);
- CommandParameter pinputdir("inputdir", "String", "", "", "", "", "",false,false); parameters.push_back(pinputdir);
- CommandParameter poutputdir("outputdir", "String", "", "", "", "", "",false,false); parameters.push_back(poutputdir);
- CommandParameter pprocessors("processors", "Number", "", "1", "", "", "",false,false); parameters.push_back(pprocessors);
+ CommandParameter piters("iters", "Number", "", "1000", "", "", "","",false,false); parameters.push_back(piters);
+ CommandParameter pdesign("design", "InputTypes", "", "", "TreeDesign", "TreeDesign", "none","summary",false,false,true); parameters.push_back(pdesign);
+ CommandParameter pshared("shared", "InputTypes", "", "", "SharedRel", "SharedRel", "none","summary",false,false,true); parameters.push_back(pshared);
+ CommandParameter prelabund("relabund", "InputTypes", "", "", "SharedRel", "SharedRel", "none","summary",false,false); parameters.push_back(prelabund);
+ CommandParameter pgroups("groups", "String", "", "", "", "", "","",false,false); parameters.push_back(pgroups);
+ CommandParameter plabel("label", "String", "", "", "", "", "","",false,false); parameters.push_back(plabel);
+ CommandParameter ptree("tree", "InputTypes", "", "", "TreeDesign", "TreeDesign", "none","tree-summary",false,false,true); parameters.push_back(ptree);
+ CommandParameter pinputdir("inputdir", "String", "", "", "", "", "","",false,false); parameters.push_back(pinputdir);
+ CommandParameter poutputdir("outputdir", "String", "", "", "", "", "","",false,false); parameters.push_back(poutputdir);
+ CommandParameter pprocessors("processors", "Number", "", "1", "", "", "","",false,false); parameters.push_back(pprocessors);
vector<string> myArray;
for (int i = 0; i < parameters.size(); i++) { myArray.push_back(parameters[i].name); }
exit(1);
}
}
-
+//**********************************************************************************************************************
+string IndicatorCommand::getOutputPattern(string type) {
+ try {
+ string pattern = "";
+
+ if (type == "tree") { pattern = "[filename],indicator.tre"; }
+ else if (type == "summary") { pattern = "[filename],indicator.summary"; }
+ else { m->mothurOut("[ERROR]: No definition for type " + type + " output pattern.\n"); m->control_pressed = true; }
+
+ return pattern;
+ }
+ catch(exception& e) {
+ m->errorOut(e, "IndicatorCommand", "getOutputPattern");
+ exit(1);
+ }
+}
//**********************************************************************************************************************
IndicatorCommand::IndicatorCommand(){
try {
m->clearGroups();
m->clearAllGroups();
m->Treenames.clear();
- m->names.clear();
vector<string> tempOutNames;
outputTypes["tree"] = tempOutNames;
designMap->readDesignMap();
//fill Groups - checks for "all" and for any typo groups
- SharedUtil* util = new SharedUtil();
+ SharedUtil util;
vector<string> nameGroups = designMap->getNamesOfGroups();
- util->setGroups(Groups, nameGroups);
+ util.setGroups(Groups, nameGroups);
designMap->setNamesOfGroups(nameGroups);
- delete util;
- //loop through the Groups and fill Globaldata's Groups with the design file info
vector<string> namesSeqs = designMap->getNamesSeqs(Groups);
m->setGroups(namesSeqs);
}
string groupfile = "";
m->setTreeFile(treefile);
Tree* tree = new Tree(treefile); delete tree; //extracts names from tree to make faked out groupmap
- treeMap = new TreeMap();
+ ct = new CountTable();
bool mismatch = false;
-
- for (int i = 0; i < m->Treenames.size(); i++) {
- //sanity check - is this a group that is not in the sharedfile?
+
+ set<string> nameMap;
+ map<string, string> groupMap;
+ set<string> gps;
+ for (int i = 0; i < m->Treenames.size(); i++) {
+ nameMap.insert(m->Treenames[i]);
+ //sanity check - is this a group that is not in the sharedfile?
+ if (i == 0) { gps.insert("Group1"); }
if (designfile == "") {
if (!(m->inUsersGroups(m->Treenames[i], m->getAllGroups()))) {
m->mothurOut("[ERROR]: " + m->Treenames[i] + " is not a group in your shared or relabund file."); m->mothurOutEndLine();
mismatch = true;
}
- treeMap->addSeq(m->Treenames[i], "Group1");
+ groupMap[m->Treenames[i]] = "Group1";
}else{
vector<string> myGroups; myGroups.push_back(m->Treenames[i]);
vector<string> myNames = designMap->getNamesSeqs(myGroups);
mismatch = true;
}
}
- treeMap->addSeq(m->Treenames[i], "Group1");
+ groupMap[m->Treenames[i]] = "Group1";
}
- }
+ }
+ ct->createTable(nameMap, groupMap, gps);
if ((designfile != "") && (m->Treenames.size() != Groups.size())) { cout << Groups.size() << '\t' << m->Treenames.size() << endl; m->mothurOut("[ERROR]: You design file does not match your tree, aborting."); m->mothurOutEndLine(); mismatch = true; }
if (designfile != "") { delete designMap; }
if (sharedfile != "") { for (int i = 0; i < lookup.size(); i++) { delete lookup[i]; } }
else { for (int i = 0; i < lookupFloat.size(); i++) { delete lookupFloat[i]; } }
- delete treeMap;
+ delete ct;
return 0;
}
read = new ReadNewickTree(treefile);
- int readOk = read->read(treeMap);
+ int readOk = read->read(ct);
- if (readOk != 0) { m->mothurOut("Read Terminated."); m->mothurOutEndLine(); delete treeMap; delete read; return 0; }
+ if (readOk != 0) { m->mothurOut("Read Terminated."); m->mothurOutEndLine(); delete ct; delete read; return 0; }
vector<Tree*> T = read->getTrees();
if (designfile != "") { delete designMap; }
if (sharedfile != "") { for (int i = 0; i < lookup.size(); i++) { delete lookup[i]; } }
else { for (int i = 0; i < lookupFloat.size(); i++) { delete lookupFloat[i]; } }
- for (int i = 0; i < T.size(); i++) { delete T[i]; } delete treeMap; return 0;
+ for (int i = 0; i < T.size(); i++) { delete T[i]; } delete ct; return 0;
}
-
+
T[0]->assembleTree();
/***************************************************/
// create ouptut tree - respecting pickedGroups //
/***************************************************/
- Tree* outputTree = new Tree(m->getNumGroups(), treeMap);
+ Tree* outputTree = new Tree(m->getNumGroups(), ct);
outputTree->getSubTree(T[0], m->getGroups());
outputTree->assembleTree();
//no longer need original tree, we have output tree to use and label
for (int i = 0; i < T.size(); i++) { delete T[i]; }
-
if (m->control_pressed) {
if (designfile != "") { delete designMap; }
if (sharedfile != "") { for (int i = 0; i < lookup.size(); i++) { delete lookup[i]; } }
else { for (int i = 0; i < lookupFloat.size(); i++) { delete lookupFloat[i]; } }
- delete outputTree; delete treeMap; return 0;
+ delete outputTree; delete ct; return 0;
}
/***************************************************/
// get indicator species values //
/***************************************************/
GetIndicatorSpecies(outputTree);
- delete outputTree; delete treeMap;
+ delete outputTree; delete ct;
}else { //run with design file only
//get indicator species
try {
string thisOutputDir = outputDir;
if (outputDir == "") { thisOutputDir += m->hasPath(inputFileName); }
- string outputFileName = thisOutputDir + m->getRootName(m->getSimpleName(inputFileName)) + "indicator.summary";
+ map<string, string> variables;
+ variables["[filename]"] = thisOutputDir + m->getRootName(m->getSimpleName(inputFileName));
+ string outputFileName = getOutputFileName("summary", variables);
outputNames.push_back(outputFileName); outputTypes["summary"].push_back(outputFileName);
ofstream out;
m->openOutputFile(outputFileName, out);
out.setf(ios::fixed, ios::floatfield); out.setf(ios::showpoint);
- m->mothurOutEndLine(); m->mothurOut("Species\tIndicatorValue\tpValue\n");
+ m->mothurOutEndLine(); m->mothurOut("Species\tIndicator_Groups\tIndicatorValue\tpValue\n");
int numBins = 0;
if (sharedfile != "") { numBins = lookup[0]->getNumBins(); }
vector<float> indicatorValues; //size of numBins
vector<float> pValues;
+ vector<string> indicatorGroups;
map< vector<int>, vector<int> > randomGroupingsMap; //maps location in groupings to location in groupings, ie, [0][0] -> [1][2]. This is so we don't have to actually move the sharedRabundVectors.
if (sharedfile != "") {
if (groupsAlreadyAdded.size() != lookup.size()) { m->mothurOut("[ERROR]: could not make proper groupings."); m->mothurOutEndLine(); }
- indicatorValues = getValues(groupings, randomGroupingsMap);
+ indicatorValues = getValues(groupings, indicatorGroups, randomGroupingsMap);
- pValues = getPValues(groupings, randomGroupingsMap, lookup.size(), indicatorValues);
+ pValues = getPValues(groupings, lookup.size(), indicatorValues);
}else {
vector< vector<SharedRAbundFloatVector*> > groupings;
set<string> groupsAlreadyAdded;
if (groupsAlreadyAdded.size() != lookupFloat.size()) { m->mothurOut("[ERROR]: could not make proper groupings."); m->mothurOutEndLine(); }
- indicatorValues = getValues(groupings, randomGroupingsMap);
+ indicatorValues = getValues(groupings, indicatorGroups, randomGroupingsMap);
- pValues = getPValues(groupings, randomGroupingsMap, lookupFloat.size(), indicatorValues);
+ pValues = getPValues(groupings, lookupFloat.size(), indicatorValues);
}
if (m->control_pressed) { out.close(); return 0; }
/******************************************************/
//output indicator values to table form //
/*****************************************************/
- out << "OTU\tIndicator_Value\tpValue" << endl;
+ out << "OTU\tIndicator_Groups\tIndicator_Value\tpValue" << endl;
for (int j = 0; j < indicatorValues.size(); j++) {
if (m->control_pressed) { out.close(); return 0; }
- out << (j+1) << '\t' << indicatorValues[j] << '\t';
+ out << m->currentBinLabels[j] << '\t' << indicatorGroups[j] << '\t' << indicatorValues[j] << '\t';
if (pValues[j] > (1/(float)iters)) { out << pValues[j] << endl; }
else { out << "<" << (1/(float)iters) << endl; }
if (pValues[j] <= 0.05) {
- cout << "OTU" << j+1 << '\t' << indicatorValues[j] << '\t';
+ cout << m->currentBinLabels[j] << '\t' << indicatorGroups[j] << '\t' << indicatorValues[j] << '\t';
string pValueString = "<" + toString((1/(float)iters));
if (pValues[j] > (1/(float)iters)) { pValueString = toString(pValues[j]); cout << pValues[j];}
else { cout << "<" << (1/(float)iters); }
- m->mothurOutJustToLog("OTU" + toString(j+1) + "\t" + toString(indicatorValues[j]) + "\t" + pValueString);
+ m->mothurOutJustToLog(m->currentBinLabels[j] + "\t" + indicatorGroups[j] + "\t" + toString(indicatorValues[j]) + "\t" + pValueString);
m->mothurOutEndLine();
}
}
string thisOutputDir = outputDir;
if (outputDir == "") { thisOutputDir += m->hasPath(inputFileName); }
- string outputFileName = thisOutputDir + m->getRootName(m->getSimpleName(inputFileName)) + "indicator.summary";
+ map<string, string> variables;
+ variables["[filename]"] = thisOutputDir + m->getRootName(m->getSimpleName(inputFileName));
+ string outputFileName = getOutputFileName("summary",variables);
outputNames.push_back(outputFileName); outputTypes["summary"].push_back(outputFileName);
ofstream out;
//print headings
out << "TreeNode\t";
- for (int i = 0; i < numBins; i++) { out << "OTU" << (i+1) << "_IndValue" << '\t' << "pValue" << '\t'; }
+ for (int i = 0; i < numBins; i++) { out << m->currentBinLabels[i] << "_IndGroups" << '\t' << m->currentBinLabels[i] << "_IndValue" << '\t' << "pValue" << '\t'; }
out << endl;
- m->mothurOutEndLine(); m->mothurOut("Node\tSpecies\tIndicatorValue\tpValue\n");
+ m->mothurOutEndLine(); m->mothurOut("Node\tSpecies\tIndicator_Groups\tIndicatorValue\tpValue\n");
string treeOutputDir = outputDir;
if (outputDir == "") { treeOutputDir += m->hasPath(treefile); }
- string outputTreeFileName = treeOutputDir + m->getRootName(m->getSimpleName(treefile)) + "indicator.tre";
+ variables["[filename]"] = treeOutputDir + m->getRootName(m->getSimpleName(treefile));
+ string outputTreeFileName = getOutputFileName("tree", variables);
//create a map from tree node index to names of descendants, save time later to know which sharedRabund you need
vector<float> indicatorValues; //size of numBins
vector<float> pValues;
+ vector<string> indicatorGroups;
map< vector<int>, vector<int> > randomGroupingsMap; //maps location in groupings to location in groupings, ie, [0][0] -> [1][2]. This is so we don't have to actually move the sharedRabundVectors.
if (sharedfile != "") {
if (groupsAlreadyAdded.size() != lookup.size()) { m->mothurOut("[ERROR]: could not make proper groupings."); m->mothurOutEndLine(); }
- indicatorValues = getValues(groupings, randomGroupingsMap);
+ indicatorValues = getValues(groupings, indicatorGroups, randomGroupingsMap);
- pValues = getPValues(groupings, randomGroupingsMap, lookup.size(), indicatorValues);
+ pValues = getPValues(groupings, lookup.size(), indicatorValues);
}else {
vector< vector<SharedRAbundFloatVector*> > groupings;
if (groupsAlreadyAdded.size() != lookupFloat.size()) { m->mothurOut("[ERROR]: could not make proper groupings."); m->mothurOutEndLine(); }
- indicatorValues = getValues(groupings, randomGroupingsMap);
+ indicatorValues = getValues(groupings, indicatorGroups, randomGroupingsMap);
- pValues = getPValues(groupings, randomGroupingsMap, lookupFloat.size(), indicatorValues);
+ pValues = getPValues(groupings, lookupFloat.size(), indicatorValues);
}
if (m->control_pressed) { out.close(); return 0; }
if (m->control_pressed) { out.close(); return 0; }
if (pValues[j] < (1/(float)iters)) {
- out << indicatorValues[j] << '\t' << '<' << (1/(float)iters) << '\t';
+ out << indicatorGroups[j] << '\t' << indicatorValues[j] << '\t' << '<' << (1/(float)iters) << '\t';
}else {
- out << indicatorValues[j] << '\t' << pValues[j] << '\t';
+ out << indicatorGroups[j] << '\t' << indicatorValues[j] << '\t' << pValues[j] << '\t';
}
if (pValues[j] <= 0.05) {
- cout << i+1 << "\tOTU" << j+1 << '\t' << indicatorValues[j] << '\t';
+ cout << i+1 << '\t' << m->currentBinLabels[j] << '\t' << indicatorGroups[j] << '\t' << indicatorValues[j] << '\t';
string pValueString = "<" + toString((1/(float)iters));
if (pValues[j] > (1/(float)iters)) { pValueString = toString(pValues[j]); cout << pValues[j];}
else { cout << "<" << (1/(float)iters); }
- m->mothurOutJustToLog(toString(i) + "\tOTU" + toString(j+1) + "\t" + toString(indicatorValues[j]) + "\t" + pValueString);
+ m->mothurOutJustToLog(toString(i) + "\t" + m->currentBinLabels[j] + "\t" + indicatorGroups[j] + "\t" + toString(indicatorValues[j]) + "\t" + pValueString);
m->mothurOutEndLine();
}
}
}
}
//**********************************************************************************************************************
-vector<float> IndicatorCommand::getValues(vector< vector<SharedRAbundFloatVector*> >& groupings, map< vector<int>, vector<int> > groupingsMap){
+vector<float> IndicatorCommand::getValues(vector< vector<SharedRAbundFloatVector*> >& groupings, vector<string>& indicatorGroupings, map< vector<int>, vector<int> > groupingsMap){
try {
vector<float> values;
map< vector<int>, vector<int> >::iterator it;
-
+
+ indicatorGroupings.clear();
+
+ //create grouping strings
+ vector<string> groupingsGroups;
+ for (int j = 0; j < groupings.size(); j++) {
+ string tempGrouping = "";
+ for (int k = 0; k < groupings[j].size()-1; k++) {
+ tempGrouping += groupings[j][k]->getGroup() + "-";
+ }
+ tempGrouping += groupings[j][groupings[j].size()-1]->getGroup();
+ groupingsGroups.push_back(tempGrouping);
+ }
+
+
//for each otu
for (int i = 0; i < groupings[0][0]->getNumBins(); i++) {
}
float maxIndVal = 0.0;
+ string maxGrouping = "";
for (int j = 0; j < terms.size(); j++) {
float thisAij = (terms[j] / AijDenominator); //relative abundance
float thisValue = thisAij * Bij[j] * 100.0;
//save largest
- if (thisValue > maxIndVal) { maxIndVal = thisValue; }
+ if (thisValue > maxIndVal) { maxIndVal = thisValue; maxGrouping = groupingsGroups[j]; }
}
values.push_back(maxIndVal);
+ indicatorGroupings.push_back(maxGrouping);
}
return values;
}
//**********************************************************************************************************************
//same as above, just data type difference
-vector<float> IndicatorCommand::getValues(vector< vector<SharedRAbundVector*> >& groupings, map< vector<int>, vector<int> > groupingsMap){
+vector<float> IndicatorCommand::getValues(vector< vector<SharedRAbundVector*> >& groupings, vector<string>& indicatorGroupings, map< vector<int>, vector<int> > groupingsMap){
try {
vector<float> values;
-
- /*for (int j = 0; j < groupings.size(); j++) {
- cout << "grouping " << j << endl;
- for (int k = 0; k < groupings[j].size(); k++) {
- cout << groupings[j][k]->getGroup() << endl;
- }
- }*/
map< vector<int>, vector<int> >::iterator it;
+
+ indicatorGroupings.clear();
+
+ //create grouping strings
+ vector<string> groupingsGroups;
+ for (int j = 0; j < groupings.size(); j++) {
+ string tempGrouping = "";
+ for (int k = 0; k < groupings[j].size()-1; k++) {
+ tempGrouping += groupings[j][k]->getGroup() + "-";
+ }
+ tempGrouping += groupings[j][groupings[j].size()-1]->getGroup();
+ groupingsGroups.push_back(tempGrouping);
+ }
+
//for each otu
for (int i = 0; i < groupings[0][0]->getNumBins(); i++) {
}
float maxIndVal = 0.0;
+ string maxGrouping = "";
for (int j = 0; j < terms.size(); j++) {
float thisAij = (terms[j] / AijDenominator); //relative abundance
float thisValue = thisAij * Bij[j] * 100.0;
//save largest
- if (thisValue > maxIndVal) { maxIndVal = thisValue; }
+ if (thisValue > maxIndVal) { maxIndVal = thisValue; maxGrouping = groupingsGroups[j]; }
}
values.push_back(maxIndVal);
+ indicatorGroupings.push_back(maxGrouping);
}
return values;
}
}
//**********************************************************************************************************************
-vector<float> IndicatorCommand::driver(vector< vector<SharedRAbundFloatVector*> >& groupings, map< vector<int>, vector<int> > groupingsMap, int num, vector<float> indicatorValues, int numIters){
+vector<float> IndicatorCommand::driver(vector< vector<SharedRAbundFloatVector*> >& groupings, int num, vector<float> indicatorValues, int numIters){
try {
vector<float> pvalues;
pvalues.resize(indicatorValues.size(), 0);
+ vector<string> notUsedGroupings; //we dont care about the grouping for the pvalues since they are randomized, but we need to pass the function something to make it work.
for(int i=0;i<numIters;i++){
if (m->control_pressed) { break; }
- groupingsMap = randomizeGroupings(groupings, num);
- vector<float> randomIndicatorValues = getValues(groupings, groupingsMap);
+ map< vector<int>, vector<int> > groupingsMap = randomizeGroupings(groupings, num);
+ vector<float> randomIndicatorValues = getValues(groupings, notUsedGroupings, groupingsMap);
for (int j = 0; j < indicatorValues.size(); j++) {
if (randomIndicatorValues[j] >= indicatorValues[j]) { pvalues[j]++; }
}
}
//**********************************************************************************************************************
-vector<float> IndicatorCommand::getPValues(vector< vector<SharedRAbundFloatVector*> >& groupings, map< vector<int>, vector<int> > groupingsMap, int num, vector<float> indicatorValues){
+vector<float> IndicatorCommand::getPValues(vector< vector<SharedRAbundFloatVector*> >& groupings, int num, vector<float> indicatorValues){
try {
vector<float> pvalues;
-
-#if defined (__APPLE__) || (__MACH__) || (linux) || (__linux)
+
if(processors == 1){
- pvalues = driver(groupings, groupingsMap, num, indicatorValues, iters);
- for (int i = 0; i < pvalues.size(); i++) { pvalues[i] /= (double)iters; }
+ pvalues = driver(groupings, num, indicatorValues, iters);
+ for (int i = 0; i < pvalues.size(); i++) { pvalues[i] /= (double)iters; }
}else{
-
- //divide iters between processors
- int numItersPerProcessor = iters / processors;
-
- vector<int> processIDS;
- int process = 1;
- int num = 0;
-
+ //divide iters between processors
+ vector<int> procIters;
+ int numItersPerProcessor = iters / processors;
+
+ //divide iters between processes
+ for (int h = 0; h < processors; h++) {
+ if(h == processors - 1){ numItersPerProcessor = iters - h * numItersPerProcessor; }
+ procIters.push_back(numItersPerProcessor);
+ }
+
+ vector<int> processIDS;
+ int process = 1;
+
+#if defined (__APPLE__) || (__MACH__) || (linux) || (__linux) || (__linux__) || (__unix__) || (__unix)
+
//loop through and create all the processes you want
while (process != processors) {
int pid = fork();
processIDS.push_back(pid); //create map from line number to pid so you can append files in correct order later
process++;
}else if (pid == 0){
- pvalues = driver(groupings, groupingsMap, num, indicatorValues, numItersPerProcessor);
+ pvalues = driver(groupings, num, indicatorValues, procIters[process]);
//pass pvalues to parent
ofstream out;
}
//do my part
- //special case for last processor in case it doesn't divide evenly
- numItersPerProcessor = iters - ((processors-1) * numItersPerProcessor);
- pvalues = driver(groupings, groupingsMap, num, indicatorValues, numItersPerProcessor);
+ pvalues = driver(groupings, num, indicatorValues, procIters[0]);
//force parent to wait until all the processes are done
for (int i=0;i<processIDS.size();i++) {
}
in.close(); m->mothurRemove(tempFile);
}
- for (int i = 0; i < pvalues.size(); i++) { pvalues[i] /= (double)iters; }
- }
+ for (int i = 0; i < pvalues.size(); i++) { pvalues[i] /= (double)iters; }
#else
- pvalues = driver(groupings, groupingsMap, num, indicatorValues, iters);
- for (int i = 0; i < pvalues.size(); i++) { pvalues[i] /= (double)iters; }
+
+ //fill in functions
+ vector<indicatorData*> pDataArray;
+ DWORD dwThreadIdArray[processors-1];
+ HANDLE hThreadArray[processors-1];
+
+ //Create processor worker threads.
+ for( int i=1; i<processors; i++ ){
+
+ //make copy of lookup so we don't get access violations
+ vector< vector<SharedRAbundFloatVector*> > newGroupings;
+
+ for (int k = 0; k < groupings.size(); k++) {
+ vector<SharedRAbundFloatVector*> newLookup;
+ for (int l = 0; l < groupings[k].size(); l++) {
+ SharedRAbundFloatVector* temp = new SharedRAbundFloatVector();
+ temp->setLabel(groupings[k][l]->getLabel());
+ temp->setGroup(groupings[k][l]->getGroup());
+ newLookup.push_back(temp);
+ }
+ newGroupings.push_back(newLookup);
+ }
+
+ //for each bin
+ for (int l = 0; l < groupings.size(); l++) {
+ for (int k = 0; k < groupings[l][0]->getNumBins(); k++) {
+ if (m->control_pressed) { for (int j = 0; j < newGroupings.size(); j++) { for (int u = 0; u < newGroupings[j].size(); u++) { delete newGroupings[j][u]; } } return pvalues; }
+
+ for (int j = 0; j < groupings[l].size(); j++) { newGroupings[l][j]->push_back(groupings[l][j]->getAbundance(k), groupings[l][j]->getGroup()); }
+ }
+ }
+
+ vector<float> copyIValues = indicatorValues;
+
+ indicatorData* temp = new indicatorData(m, procIters[i], newGroupings, num, copyIValues);
+ pDataArray.push_back(temp);
+ processIDS.push_back(i);
+
+ hThreadArray[i-1] = CreateThread(NULL, 0, MyIndicatorThreadFunction, pDataArray[i-1], 0, &dwThreadIdArray[i-1]);
+ }
+
+ //do my part
+ pvalues = driver(groupings, num, indicatorValues, procIters[0]);
+
+ //Wait until all threads have terminated.
+ WaitForMultipleObjects(processors-1, hThreadArray, TRUE, INFINITE);
+
+ //Close all thread handles and free memory allocations.
+ for(int i=0; i < pDataArray.size(); i++){
+ for (int j = 0; j < pDataArray[i]->pvalues.size(); j++) { pvalues[j] += pDataArray[i]->pvalues[j]; }
+
+ for (int l = 0; l < pDataArray[i]->groupings.size(); l++) {
+ for (int j = 0; j < pDataArray[i]->groupings[l].size(); j++) { delete pDataArray[i]->groupings[l][j]; }
+ }
+
+ CloseHandle(hThreadArray[i]);
+ delete pDataArray[i];
+ }
+
+ for (int i = 0; i < pvalues.size(); i++) { pvalues[i] /= (double)iters; }
#endif
+ }
+
return pvalues;
}
//**********************************************************************************************************************
//same as above, just data type difference
-vector<float> IndicatorCommand::driver(vector< vector<SharedRAbundVector*> >& groupings, map< vector<int>, vector<int> > groupingsMap, int num, vector<float> indicatorValues, int numIters){
+vector<float> IndicatorCommand::driver(vector< vector<SharedRAbundVector*> >& groupings, int num, vector<float> indicatorValues, int numIters){
try {
vector<float> pvalues;
pvalues.resize(indicatorValues.size(), 0);
+ vector<string> notUsedGroupings; //we dont care about the grouping for the pvalues since they are randomized, but we need to pass the function something to make it work.
for(int i=0;i<numIters;i++){
if (m->control_pressed) { break; }
- groupingsMap = randomizeGroupings(groupings, num);
- vector<float> randomIndicatorValues = getValues(groupings, groupingsMap);
+ map< vector<int>, vector<int> > groupingsMap = randomizeGroupings(groupings, num);
+ vector<float> randomIndicatorValues = getValues(groupings, notUsedGroupings, groupingsMap);
for (int j = 0; j < indicatorValues.size(); j++) {
if (randomIndicatorValues[j] >= indicatorValues[j]) { pvalues[j]++; }
}
//**********************************************************************************************************************
//same as above, just data type difference
-vector<float> IndicatorCommand::getPValues(vector< vector<SharedRAbundVector*> >& groupings, map< vector<int>, vector<int> > groupingsMap, int num, vector<float> indicatorValues){
+vector<float> IndicatorCommand::getPValues(vector< vector<SharedRAbundVector*> >& groupings, int num, vector<float> indicatorValues){
try {
vector<float> pvalues;
-
-#if defined (__APPLE__) || (__MACH__) || (linux) || (__linux)
+
if(processors == 1){
- pvalues = driver(groupings, groupingsMap, num, indicatorValues, iters);
- for (int i = 0; i < pvalues.size(); i++) { pvalues[i] /= (double)iters; }
+ pvalues = driver(groupings, num, indicatorValues, iters);
+ for (int i = 0; i < pvalues.size(); i++) { pvalues[i] /= (double)iters; }
}else{
-
- //divide iters between processors
- int numItersPerProcessor = iters / processors;
-
- vector<int> processIDS;
- int process = 1;
-
+ //divide iters between processors
+ vector<int> procIters;
+ int numItersPerProcessor = iters / processors;
+
+ //divide iters between processes
+ for (int h = 0; h < processors; h++) {
+ if(h == processors - 1){ numItersPerProcessor = iters - h * numItersPerProcessor; }
+ procIters.push_back(numItersPerProcessor);
+ }
+
+ vector<int> processIDS;
+ int process = 1;
+
+#if defined (__APPLE__) || (__MACH__) || (linux) || (__linux) || (__linux__) || (__unix__) || (__unix)
+
//loop through and create all the processes you want
while (process != processors) {
int pid = fork();
processIDS.push_back(pid); //create map from line number to pid so you can append files in correct order later
process++;
}else if (pid == 0){
- pvalues = driver(groupings, groupingsMap, num, indicatorValues, numItersPerProcessor);
+ pvalues = driver(groupings, num, indicatorValues, procIters[process]);
//pass pvalues to parent
ofstream out;
out.close();
exit(0);
- }else {
- m->mothurOut("[ERROR]: unable to spawn the necessary processes."); m->mothurOutEndLine();
+ }else {
+ m->mothurOut("[ERROR]: unable to spawn the necessary processes."); m->mothurOutEndLine();
for (int i = 0; i < processIDS.size(); i++) { kill (processIDS[i], SIGINT); }
exit(0);
}
}
//do my part
- //special case for last processor in case it doesn't divide evenly
- numItersPerProcessor = iters - ((processors-1) * numItersPerProcessor);
- pvalues = driver(groupings, groupingsMap, num, indicatorValues, numItersPerProcessor);
+ pvalues = driver(groupings, num, indicatorValues, procIters[0]);
//force parent to wait until all the processes are done
- for (int i=0;i<processIDS.size();i++) {
+ for (int i=0;i<processIDS.size();i++) {
int temp = processIDS[i];
wait(&temp);
}
- //combine results
+ //combine results
for (int i = 0; i < processIDS.size(); i++) {
ifstream in;
string tempFile = toString(processIDS[i]) + ".pvalues.temp";
m->openInputFile(tempFile, in);
////// to do ///////////
- int numTemp; numTemp = 0;
+ int numTemp; numTemp = 0;
for (int j = 0; j < pvalues.size(); j++) {
in >> numTemp; m->gobble(in);
pvalues[j] += numTemp;
}
in.close(); m->mothurRemove(tempFile);
}
- for (int i = 0; i < pvalues.size(); i++) { pvalues[i] /= (double)iters; }
- }
+ for (int i = 0; i < pvalues.size(); i++) { pvalues[i] /= (double)iters; }
#else
- pvalues = driver(groupings, groupingsMap, num, indicatorValues, iters);
- for (int i = 0; i < pvalues.size(); i++) { pvalues[i] /= (double)iters; }
+
+ //fill in functions
+ vector<indicatorData*> pDataArray;
+ DWORD dwThreadIdArray[processors-1];
+ HANDLE hThreadArray[processors-1];
+
+ //Create processor worker threads.
+ for( int i=1; i<processors; i++ ){
+
+ //make copy of lookup so we don't get access violations
+ vector< vector<SharedRAbundFloatVector*> > newGroupings;
+
+ for (int k = 0; k < groupings.size(); k++) {
+ vector<SharedRAbundFloatVector*> newLookup;
+ for (int l = 0; l < groupings[k].size(); l++) {
+ SharedRAbundFloatVector* temp = new SharedRAbundFloatVector();
+ temp->setLabel(groupings[k][l]->getLabel());
+ temp->setGroup(groupings[k][l]->getGroup());
+ newLookup.push_back(temp);
+ }
+ newGroupings.push_back(newLookup);
+ }
+
+ //for each bin
+ for (int l = 0; l < groupings.size(); l++) {
+ for (int k = 0; k < groupings[l][0]->getNumBins(); k++) {
+ if (m->control_pressed) { for (int j = 0; j < newGroupings.size(); j++) { for (int u = 0; u < newGroupings[j].size(); u++) { delete newGroupings[j][u]; } } return pvalues; }
+
+ for (int j = 0; j < groupings[l].size(); j++) { newGroupings[l][j]->push_back((float)(groupings[l][j]->getAbundance(k)), groupings[l][j]->getGroup()); }
+ }
+ }
+
+ vector<float> copyIValues = indicatorValues;
+
+ indicatorData* temp = new indicatorData(m, procIters[i], newGroupings, num, copyIValues);
+ pDataArray.push_back(temp);
+ processIDS.push_back(i);
+
+ hThreadArray[i-1] = CreateThread(NULL, 0, MyIndicatorThreadFunction, pDataArray[i-1], 0, &dwThreadIdArray[i-1]);
+ }
+
+ //do my part
+ pvalues = driver(groupings, num, indicatorValues, procIters[0]);
+
+ //Wait until all threads have terminated.
+ WaitForMultipleObjects(processors-1, hThreadArray, TRUE, INFINITE);
+
+ //Close all thread handles and free memory allocations.
+ for(int i=0; i < pDataArray.size(); i++){
+ for (int j = 0; j < pDataArray[i]->pvalues.size(); j++) { pvalues[j] += pDataArray[i]->pvalues[j]; }
+
+ for (int l = 0; l < pDataArray[i]->groupings.size(); l++) {
+ for (int j = 0; j < pDataArray[i]->groupings[l].size(); j++) { delete pDataArray[i]->groupings[l][j]; }
+ }
+
+ CloseHandle(hThreadArray[i]);
+ delete pDataArray[i];
+ }
+
+ for (int i = 0; i < pvalues.size(); i++) { pvalues[i] /= (double)iters; }
#endif
+ }
+
return pvalues;
}
catch(exception& e) {
- m->errorOut(e, "IndicatorCommand", "getPValues");
+ m->errorOut(e, "IndicatorCommand", "getPValues");
exit(1);
}
}
projects / mothur.git / blobdiff