2 // getmetacommunitycommand.h
5 // Created by SarahsWork on 4/9/13.
6 // Copyright (c) 2013 Schloss Lab. All rights reserved.
9 #ifndef Mothur_getmetacommunitycommand_h
10 #define Mothur_getmetacommunitycommand_h
12 #include "command.hpp"
13 #include "inputdata.h"
14 #include "qFinderDMM.h"
16 /**************************************************************************************************/
18 class GetMetaCommunityCommand : public Command {
20 GetMetaCommunityCommand(string);
21 GetMetaCommunityCommand();
22 ~GetMetaCommunityCommand(){}
24 vector<string> setParameters();
25 string getCommandName() { return "get.metacommunity"; }
26 string getCommandCategory() { return "OTU-Based Approaches"; }
28 string getOutputPattern(string);
30 string getHelpString();
31 string getCitation() { return "http://www.mothur.org/wiki/Get.metacommunity"; }
32 string getDescription() { return "brief description"; }
35 void help() { m->mothurOut(getHelpString()); }
39 unsigned long long start;
40 unsigned long long end;
41 linePair(unsigned long long i, unsigned long long j) : start(i), end(j) {}
45 vector<string> outputNames;
47 int minpartitions, maxpartitions, optimizegap, processors;
48 vector<string> Groups;
51 int processDriver(vector<SharedRAbundVector*>&, vector<int>&, string, vector<string>, vector<string>, vector<string>, int);
52 int createProcesses(vector<SharedRAbundVector*>&);
53 vector<double> generateDesignFile(int, map<string,string>);
54 int generateSummaryFile(int, map<string,string>);
58 /**************************************************************************************************/
62 double refMean, difference;
63 vector<double> partMean, partLCI, partUCI;
66 /**************************************************************************************************/
68 struct metaCommunityData {
69 vector<SharedRAbundVector*> thislookup;
71 string outputFileName;
72 vector<string> relabunds, matrix, outputNames;
73 int minpartitions, maxpartitions, optimizegap;
78 metaCommunityData(vector<SharedRAbundVector*> lu, MothurOut* mout, vector<int> dp, string fit, vector<string> rels, vector<string> mat, int minp, int maxp, int opg) {
91 /**************************************************************************************************/
92 #if defined (__APPLE__) || (__MACH__) || (linux) || (__linux) || (__linux__) || (__unix__) || (__unix)
94 static DWORD WINAPI MyMetaCommunityThreadFunction(LPVOID lpParam){
95 metaCommunityData* pDataArray;
96 pDataArray = (metaCommunityData*)lpParam;
100 double minLaplace = 1e10;
103 pDataArray->m->openOutputFile(pDataArray->outputFileName, fitData);
104 fitData.setf(ios::fixed, ios::floatfield);
105 fitData.setf(ios::showpoint);
106 cout.setf(ios::fixed, ios::floatfield);
107 cout.setf(ios::showpoint);
109 vector< vector<int> > sharedMatrix;
110 for (int i = 0; i < pDataArray->thislookup.size(); i++) { sharedMatrix.push_back(pDataArray->thislookup[i]->getAbundances()); }
112 pDataArray->m->mothurOut("K\tNLE\t\tlogDet\tBIC\t\tAIC\t\tLaplace\n");
113 fitData << "K\tNLE\tlogDet\tBIC\tAIC\tLaplace" << endl;
115 for(int i=0;i<pDataArray->parts.size();i++){
117 int numPartitions = pDataArray->parts[i];
119 if (pDataArray->m->debug) { pDataArray->m->mothurOut("[DEBUG]: running partition " + toString(numPartitions) + "\n"); }
121 if (pDataArray->m->control_pressed) { break; }
123 qFinderDMM* findQ = new qFinderDMM(sharedMatrix, numPartitions);
125 if (pDataArray->m->debug) { pDataArray->m->mothurOut("[DEBUG]: done finding Q " + toString(numPartitions) + "\n"); }
127 double laplace = findQ->getLaplace();
128 pDataArray->m->mothurOut(toString(numPartitions) + '\t');
129 cout << setprecision (2) << findQ->getNLL() << '\t' << findQ->getLogDet() << '\t';
130 pDataArray->m->mothurOutJustToLog(toString(findQ->getNLL()) + '\t' + toString(findQ->getLogDet()) + '\t');
131 cout << findQ->getBIC() << '\t' << findQ->getAIC() << '\t' << laplace;
132 pDataArray->m->mothurOutJustToLog(toString(findQ->getBIC()) + '\t' + toString(findQ->getAIC()) + '\t' + toString(laplace));
134 fitData << numPartitions << '\t';
135 fitData << setprecision (2) << findQ->getNLL() << '\t' << findQ->getLogDet() << '\t';
136 fitData << findQ->getBIC() << '\t' << findQ->getAIC() << '\t' << laplace << endl;
138 if(laplace < minLaplace){
139 pDataArray->minPartition = numPartitions;
140 minLaplace = laplace;
141 pDataArray->m->mothurOut("***");
143 pDataArray->m->mothurOutEndLine();
145 pDataArray->outputNames.push_back(pDataArray->relabunds[i]);
146 pDataArray->outputNames.push_back(pDataArray->matrix[i]);
148 findQ->printZMatrix(pDataArray->matrix[i], pDataArray->m->getGroups());
149 findQ->printRelAbund(pDataArray->relabunds[i], pDataArray->m->currentBinLabels);
151 if(pDataArray->optimizegap != -1 && (numPartitions - pDataArray->minPartition) >= pDataArray->optimizegap && numPartitions >= pDataArray->minpartitions){ break; }
159 if (pDataArray->m->control_pressed) { return 0; }
164 catch(exception& e) {
165 pDataArray->m->errorOut(e, "GetMetaCommunityCommand", "MyMetaCommunityThreadFunction");