372E12700F26365B0095CF7E /* readotucommand.cpp in Sources */ = {isa = PBXBuildFile; fileRef = 372E126F0F26365B0095CF7E /* readotucommand.cpp */; };
372E12960F263D5A0095CF7E /* readdistcommand.cpp in Sources */ = {isa = PBXBuildFile; fileRef = 372E12950F263D5A0095CF7E /* readdistcommand.cpp */; };
372E12ED0F264D320095CF7E /* commandfactory.cpp in Sources */ = {isa = PBXBuildFile; fileRef = 372E12EC0F264D320095CF7E /* commandfactory.cpp */; };
- 3731C1F30F8CFC0A0065A2AD /* treegroupscommand.cpp in Sources */ = {isa = PBXBuildFile; fileRef = 3731C1F20F8CFC0A0065A2AD /* treegroupscommand.cpp */; };
374610780F40645300460C57 /* unifracweightedcommand.cpp in Sources */ = {isa = PBXBuildFile; fileRef = 374610770F40645300460C57 /* unifracweightedcommand.cpp */; };
3746107E0F4064D100460C57 /* weighted.cpp in Sources */ = {isa = PBXBuildFile; fileRef = 3746107D0F4064D100460C57 /* weighted.cpp */; };
374610830F40652400460C57 /* unweighted.cpp in Sources */ = {isa = PBXBuildFile; fileRef = 374610820F40652400460C57 /* unweighted.cpp */; };
37519A6B0F80E6EB00FED5E8 /* sharedanderbergs.cpp in Sources */ = {isa = PBXBuildFile; fileRef = 37519A6A0F80E6EB00FED5E8 /* sharedanderbergs.cpp */; };
37519AA10F810D0200FED5E8 /* venncommand.cpp in Sources */ = {isa = PBXBuildFile; fileRef = 37519AA00F810D0200FED5E8 /* venncommand.cpp */; };
37519AB50F810FAE00FED5E8 /* venn.cpp in Sources */ = {isa = PBXBuildFile; fileRef = 37519AB40F810FAE00FED5E8 /* venn.cpp */; };
- 375873E70F7D63E90040F377 /* coverage.cpp in Sources */ = {isa = PBXBuildFile; fileRef = 375873E50F7D63E90040F377 /* coverage.cpp */; };
375873EC0F7D64520040F377 /* fullmatrix.cpp in Sources */ = {isa = PBXBuildFile; fileRef = 375873EB0F7D64520040F377 /* fullmatrix.cpp */; };
375873EF0F7D646F0040F377 /* heatmap.cpp in Sources */ = {isa = PBXBuildFile; fileRef = 375873EE0F7D646F0040F377 /* heatmap.cpp */; };
375873F20F7D64800040F377 /* heatmapcommand.cpp in Sources */ = {isa = PBXBuildFile; fileRef = 375873F00F7D64800040F377 /* heatmapcommand.cpp */; };
37AFC71F0F445386005F492D /* sharedsobscollectsummary.cpp in Sources */ = {isa = PBXBuildFile; fileRef = 37AFC71E0F445386005F492D /* sharedsobscollectsummary.cpp */; };
37B28F680F27590100808A62 /* deconvolutecommand.cpp in Sources */ = {isa = PBXBuildFile; fileRef = 37B28F670F27590100808A62 /* deconvolutecommand.cpp */; };
37C1D9730F86506E0059E3F0 /* binsequencecommand.cpp in Sources */ = {isa = PBXBuildFile; fileRef = 37C1D9720F86506E0059E3F0 /* binsequencecommand.cpp */; };
- 37C7F3A90F8E90AD00E91C2B /* sharedutilities.cpp in Sources */ = {isa = PBXBuildFile; fileRef = 37C7F3A80F8E90AD00E91C2B /* sharedutilities.cpp */; };
37D928550F21331F001D4494 /* ace.cpp in Sources */ = {isa = PBXBuildFile; fileRef = 37D927B80F21331F001D4494 /* ace.cpp */; };
37D928560F21331F001D4494 /* averagelinkage.cpp in Sources */ = {isa = PBXBuildFile; fileRef = 37D927BA0F21331F001D4494 /* averagelinkage.cpp */; };
37D928570F21331F001D4494 /* bootstrap.cpp in Sources */ = {isa = PBXBuildFile; fileRef = 37D927BB0F21331F001D4494 /* bootstrap.cpp */; };
37D9289F0F21331F001D4494 /* validparameter.cpp in Sources */ = {isa = PBXBuildFile; fileRef = 37D928530F21331F001D4494 /* validparameter.cpp */; };
37E5F3E30F29FD4200F8D827 /* treenode.cpp in Sources */ = {isa = PBXBuildFile; fileRef = 37E5F3E20F29FD4200F8D827 /* treenode.cpp */; };
37E5F4920F2A3DA800F8D827 /* readtreecommand.cpp in Sources */ = {isa = PBXBuildFile; fileRef = 37E5F4910F2A3DA800F8D827 /* readtreecommand.cpp */; };
+ 7E412F490F8D21B600381DD0 /* slibshuff.cpp in Sources */ = {isa = PBXBuildFile; fileRef = 7E412F480F8D21B600381DD0 /* slibshuff.cpp */; };
+ 7E412FEA0F8D3E2C00381DD0 /* libshuff.cpp in Sources */ = {isa = PBXBuildFile; fileRef = 7E412FE90F8D3E2C00381DD0 /* libshuff.cpp */; };
+ 7E4130F80F8E58FA00381DD0 /* dlibshuff.cpp in Sources */ = {isa = PBXBuildFile; fileRef = 7E4130F60F8E58FA00381DD0 /* dlibshuff.cpp */; };
8DD76F6A0486A84900D96B5E /* Mothur.1 in CopyFiles */ = {isa = PBXBuildFile; fileRef = C6859E8B029090EE04C91782 /* Mothur.1 */; };
A70B53AA0F4CD7AD0064797E /* getgroupcommand.cpp in Sources */ = {isa = PBXBuildFile; fileRef = A70B53A40F4CD7AD0064797E /* getgroupcommand.cpp */; };
A70B53AB0F4CD7AD0064797E /* getlabelcommand.cpp in Sources */ = {isa = PBXBuildFile; fileRef = A70B53A60F4CD7AD0064797E /* getlabelcommand.cpp */; };
372E12940F263D5A0095CF7E /* readdistcommand.h */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.c.h; path = readdistcommand.h; sourceTree = "<group>"; };
372E12950F263D5A0095CF7E /* readdistcommand.cpp */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.cpp.cpp; path = readdistcommand.cpp; sourceTree = "<group>"; };
372E12EC0F264D320095CF7E /* commandfactory.cpp */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.cpp.cpp; path = commandfactory.cpp; sourceTree = "<group>"; };
- 3731C1F10F8CFC0A0065A2AD /* treegroupscommand.h */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.c.h; path = treegroupscommand.h; sourceTree = "<group>"; };
- 3731C1F20F8CFC0A0065A2AD /* treegroupscommand.cpp */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.cpp.cpp; path = treegroupscommand.cpp; sourceTree = "<group>"; };
374610760F40645300460C57 /* unifracweightedcommand.h */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.c.h; path = unifracweightedcommand.h; sourceTree = "<group>"; };
374610770F40645300460C57 /* unifracweightedcommand.cpp */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.cpp.cpp; path = unifracweightedcommand.cpp; sourceTree = "<group>"; };
3746107C0F4064D100460C57 /* weighted.h */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.c.h; path = weighted.h; sourceTree = "<group>"; };
37519AA00F810D0200FED5E8 /* venncommand.cpp */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.cpp.cpp; path = venncommand.cpp; sourceTree = "<group>"; };
37519AB30F810FAE00FED5E8 /* venn.h */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.c.h; path = venn.h; sourceTree = "<group>"; };
37519AB40F810FAE00FED5E8 /* venn.cpp */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.cpp.cpp; path = venn.cpp; sourceTree = "<group>"; };
- 375873E50F7D63E90040F377 /* coverage.cpp */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.cpp.cpp; path = coverage.cpp; sourceTree = "<group>"; };
- 375873E60F7D63E90040F377 /* coverage.h */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.c.h; path = coverage.h; sourceTree = "<group>"; };
375873EA0F7D64520040F377 /* fullmatrix.h */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.c.h; path = fullmatrix.h; sourceTree = "<group>"; };
375873EB0F7D64520040F377 /* fullmatrix.cpp */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.cpp.cpp; path = fullmatrix.cpp; sourceTree = "<group>"; };
375873ED0F7D646F0040F377 /* heatmap.h */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.c.h; path = heatmap.h; sourceTree = "<group>"; };
37B28F670F27590100808A62 /* deconvolutecommand.cpp */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.cpp.cpp; path = deconvolutecommand.cpp; sourceTree = "<group>"; };
37C1D9710F86506E0059E3F0 /* binsequencecommand.h */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.c.h; path = binsequencecommand.h; sourceTree = "<group>"; };
37C1D9720F86506E0059E3F0 /* binsequencecommand.cpp */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.cpp.cpp; path = binsequencecommand.cpp; sourceTree = "<group>"; };
- 37C7F3A70F8E90AD00E91C2B /* sharedutilities.h */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.c.h; path = sharedutilities.h; sourceTree = "<group>"; };
- 37C7F3A80F8E90AD00E91C2B /* sharedutilities.cpp */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.cpp.cpp; path = sharedutilities.cpp; sourceTree = "<group>"; };
37D927B80F21331F001D4494 /* ace.cpp */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.cpp.cpp; path = ace.cpp; sourceTree = "<group>"; };
37D927B90F21331F001D4494 /* ace.h */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.c.h; path = ace.h; sourceTree = "<group>"; };
37D927BA0F21331F001D4494 /* averagelinkage.cpp */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.cpp.cpp; path = averagelinkage.cpp; sourceTree = "<group>"; };
37D928490F21331F001D4494 /* summarydisplay.h */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.c.h; path = summarydisplay.h; sourceTree = "<group>"; };
37D9284A0F21331F001D4494 /* summarysharedcommand.cpp */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.cpp.cpp; path = summarysharedcommand.cpp; sourceTree = "<group>"; };
37D9284B0F21331F001D4494 /* summarysharedcommand.h */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.c.h; path = summarysharedcommand.h; sourceTree = "<group>"; };
+ 37D9284C0F21331F001D4494 /* utilities.hpp */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.cpp.h; path = utilities.hpp; sourceTree = "<group>"; };
37D9284D0F21331F001D4494 /* uvest.cpp */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.cpp.cpp; path = uvest.cpp; sourceTree = "<group>"; };
37D9284E0F21331F001D4494 /* uvest.h */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.c.h; path = uvest.h; sourceTree = "<group>"; };
37D9284F0F21331F001D4494 /* validcalculator.cpp */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.cpp.cpp; path = validcalculator.cpp; sourceTree = "<group>"; };
37E5F3E20F29FD4200F8D827 /* treenode.cpp */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.cpp.cpp; path = treenode.cpp; sourceTree = "<group>"; };
37E5F4900F2A3DA800F8D827 /* readtreecommand.h */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.c.h; path = readtreecommand.h; sourceTree = "<group>"; };
37E5F4910F2A3DA800F8D827 /* readtreecommand.cpp */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.cpp.cpp; path = readtreecommand.cpp; sourceTree = "<group>"; };
+ 7E412F420F8D213C00381DD0 /* libshuff.h */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.c.h; path = libshuff.h; sourceTree = "<group>"; };
+ 7E412F470F8D21B600381DD0 /* slibshuff.h */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.c.h; path = slibshuff.h; sourceTree = "<group>"; };
+ 7E412F480F8D21B600381DD0 /* slibshuff.cpp */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.cpp.cpp; path = slibshuff.cpp; sourceTree = "<group>"; };
+ 7E412FE90F8D3E2C00381DD0 /* libshuff.cpp */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.cpp.cpp; path = libshuff.cpp; sourceTree = "<group>"; };
+ 7E4130F60F8E58FA00381DD0 /* dlibshuff.cpp */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.cpp.cpp; path = dlibshuff.cpp; sourceTree = "<group>"; };
+ 7E4130F70F8E58FA00381DD0 /* dlibshuff.h */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.c.h; path = dlibshuff.h; sourceTree = "<group>"; };
8DD76F6C0486A84900D96B5E /* mothur */ = {isa = PBXFileReference; explicitFileType = "compiled.mach-o.executable"; includeInIndex = 0; path = mothur; sourceTree = BUILT_PRODUCTS_DIR; };
A70B53A40F4CD7AD0064797E /* getgroupcommand.cpp */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.cpp.cpp; path = getgroupcommand.cpp; sourceTree = "<group>"; };
A70B53A50F4CD7AD0064797E /* getgroupcommand.h */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.c.h; path = getgroupcommand.h; sourceTree = "<group>"; };
08FB7795FE84155DC02AAC07 /* Source */ = {
isa = PBXGroup;
children = (
+ 7E4130F60F8E58FA00381DD0 /* dlibshuff.cpp */,
+ 7E4130F70F8E58FA00381DD0 /* dlibshuff.h */,
37D927BA0F21331F001D4494 /* averagelinkage.cpp */,
37D928A60F2133C0001D4494 /* calculators */,
37D927C20F21331F001D4494 /* cluster.hpp */,
37D927C10F21331F001D4494 /* cluster.cpp */,
+ 7E412F470F8D21B600381DD0 /* slibshuff.h */,
+ 7E412F480F8D21B600381DD0 /* slibshuff.cpp */,
+ 7E412F420F8D213C00381DD0 /* libshuff.h */,
+ 7E412FE90F8D3E2C00381DD0 /* libshuff.cpp */,
37D928A90F2133E5001D4494 /* commands */,
37D927C60F21331F001D4494 /* collect.h */,
37D927C50F21331F001D4494 /* collect.cpp */,
37D9281C0F21331F001D4494 /* sequence.cpp */,
37D928210F21331F001D4494 /* shared.h */,
37D928200F21331F001D4494 /* shared.cpp */,
- 37C7F3A70F8E90AD00E91C2B /* sharedutilities.h */,
- 37C7F3A80F8E90AD00E91C2B /* sharedutilities.cpp */,
37D928420F21331F001D4494 /* singlelinkage.cpp */,
37D928480F21331F001D4494 /* summarydata.h */,
37D928490F21331F001D4494 /* summarydisplay.h */,
+ 37D9284C0F21331F001D4494 /* utilities.hpp */,
37519AB30F810FAE00FED5E8 /* venn.h */,
37519AB40F810FAE00FED5E8 /* venn.cpp */,
);
EB1216E40F61ACFB004A865F /* bstick.cpp */,
37D927C00F21331F001D4494 /* chao1.h */,
37D927BF0F21331F001D4494 /* chao1.cpp */,
- 375873E60F7D63E90040F377 /* coverage.h */,
- 375873E50F7D63E90040F377 /* coverage.cpp */,
EB9303F70F53517300E8EF26 /* geom.h */,
EB9303F80F53517300E8EF26 /* geom.cpp */,
EB9303E90F534D9400E8EF26 /* logsd.h */,
37D928460F21331F001D4494 /* summarycommand.cpp */,
37D9284B0F21331F001D4494 /* summarysharedcommand.h */,
37D9284A0F21331F001D4494 /* summarysharedcommand.cpp */,
- 3731C1F10F8CFC0A0065A2AD /* treegroupscommand.h */,
- 3731C1F20F8CFC0A0065A2AD /* treegroupscommand.cpp */,
3746109B0F40657600460C57 /* unifracunweightedcommand.h */,
3746109C0F40657600460C57 /* unifracunweightedcommand.cpp */,
374610760F40645300460C57 /* unifracweightedcommand.h */,
EB1216E50F61ACFB004A865F /* bstick.cpp in Sources */,
EB1217230F61C9AC004A865F /* sharedkstest.cpp in Sources */,
EB6F015B0F6AC1670048081A /* sharedbdiversity.cpp in Sources */,
- 375873E70F7D63E90040F377 /* coverage.cpp in Sources */,
375873EC0F7D64520040F377 /* fullmatrix.cpp in Sources */,
375873EF0F7D646F0040F377 /* heatmap.cpp in Sources */,
375873F20F7D64800040F377 /* heatmapcommand.cpp in Sources */,
37519AB50F810FAE00FED5E8 /* venn.cpp in Sources */,
37C1D9730F86506E0059E3F0 /* binsequencecommand.cpp in Sources */,
370B88070F8A4EE4005AB382 /* getoturepcommand.cpp in Sources */,
- 3731C1F30F8CFC0A0065A2AD /* treegroupscommand.cpp in Sources */,
- 37C7F3A90F8E90AD00E91C2B /* sharedutilities.cpp in Sources */,
+ 7E412F490F8D21B600381DD0 /* slibshuff.cpp in Sources */,
+ 7E412FEA0F8D3E2C00381DD0 /* libshuff.cpp in Sources */,
+ 7E4130F80F8E58FA00381DD0 /* dlibshuff.cpp in Sources */,
);
runOnlyForDeploymentPostprocessing = 0;
};
ppc,
i386,
);
- GCC_OPTIMIZATION_LEVEL = 0;
+ GCC_OPTIMIZATION_LEVEL = 3;
GCC_WARN_ABOUT_RETURN_TYPE = YES;
GCC_WARN_UNUSED_VARIABLE = YES;
PREBINDING = NO;
+++ /dev/null
-/*
- * coverage.cpp
- * Mothur
- *
- * Created by Sarah Westcott on 3/9/09.
- * Copyright 2009 Schloss Lab UMASS Amherst. All rights reserved.
- *
- */
-
-/* This class library coverage at the given distances of the Cramer-von Mises statistic.
- you may refer to the "Integration of Microbial Ecology and Statistics: A Test To Compare Gene Libraries"
- paper in Applied and Environmental Microbiology, Sept. 2004, p. 5485-5492 0099-2240/04/$8.00+0
- DOI: 10.1128/AEM.70.9.5485-5492.2004 Copyright 2004 American Society for Microbiology for more information. */
-
-
-#include "coverage.h"
-
-//**********************************************************************************************************************
-Coverage::Coverage() {
- globaldata = GlobalData::getInstance();
- numUserGroups = globaldata->Groups.size();
- numGroups = globaldata->gGroupmap->getNumGroups();
-}
-
-//**********************************************************************************************************************
-void Coverage::getValues(FullMatrix* matrix, vector< vector< vector<float> > >& data, vector<float> dist) {
- try {
- vector<float> min;
- vector<string> groups;
-
- //initialize data
- data.resize(dist.size());
- for (int l = 0; l < data.size(); l++) {
- data[l].resize(numGroups);
- for (int k = 0; k < data[l].size(); k++) {
- data[l][k].push_back(0.0);
- }
- }
-
- /**************************************/
- //get the minimums for each comparision
- /**************************************/
- int count = 0;
-
- for (int i = 0; i < numGroups; i++) {
- for (int j = 0; j < numGroups; j++) {
-
- //is this "box" one the user wants analyzed?
- if ((inUsersGroups(globaldata->gGroupmap->namesOfGroups[i], globaldata->Groups) == true) && (inUsersGroups(globaldata->gGroupmap->namesOfGroups[j], globaldata->Groups) == true)) {
-
- min = matrix->getMins(count); //returns vector of mins for "box" requested ie. groups A, B, 0 = AA, 1 = AB, 2 = BA, 3 = BB;
-
- //find the coverage at this distance
- sort(min.begin(), min.end());
-
- //loop through each distance and fill data
- for (int k = 0; k < data.size(); k++) {
-
- int index = -1;
- //find index in min where value is higher than d
- for (int m = 0; m < min.size(); m++) {
- if (min[m] > dist[k]) { index = m; break; }
- }
-
- // if you don't find one than all the mins are less than d
- if (index == -1) { index = min.size(); }
-
- //save value in data
- data[k][i][j] = 1.0 - ((min.size()-index)/(float)min.size());
-
- }
- }
- count++;
- }
- }
-
- }
- catch(exception& e) {
- cout << "Standard Error: " << e.what() << " has occurred in the Coverage class Function getValues. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
- exit(1);
- }
- catch(...) {
- cout << "An unknown error has occurred in the Coverage class function getValues. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
- exit(1);
- }
-}
-
-//**********************************************************************************************************************
-void Coverage::getValues(FullMatrix* matrix, vector< vector< vector<float> > >& data, vector<float> dist, string mode) {
- try {
- vector<float> min1;
- vector<float> min2;
- vector<string> groups;
-
- //initialize data
- data.resize(dist.size());
- for (int l = 0; l < data.size(); l++) {
- data[l].resize(numGroups);
- for (int k = 0; k < data[l].size(); k++) {
- data[l][k].push_back(0.0);
- }
- }
-
- int count = 0;
- int count2 = 0;
-
- //for each box
- for (int i = 0; i < numGroups; i++) {
- for (int j = 0; j < numGroups; j++) {
-
- if (i != j) {
- //is this "box" one the user wants analyzed?
- if ((inUsersGroups(globaldata->gGroupmap->namesOfGroups[i], globaldata->Groups) == true) && (inUsersGroups(globaldata->gGroupmap->namesOfGroups[j], globaldata->Groups) == true)) {
-
- matrix->shuffle(globaldata->gGroupmap->namesOfGroups[i], globaldata->gGroupmap->namesOfGroups[j]);
-
- min1 = matrix->getMins(count); //returns vector of mins for "box" requested ie. groups A, B, 0 = AA, 1 = AB, 2 = BA, 3 = BB;
- min2 = matrix->getMins(count2); //returns vector of mins for "box" requested ie. groups A, B, 0 = AA, 1 = AB, 2 = BA, 3 = BB;
-
- //find the coverage at this distance
- sort(min1.begin(), min1.end());
-
- //find the coverage at this distance
- sort(min2.begin(), min2.end());
-
- float distDiff = 0;
-
- //loop through each distance and fill data
- for (int k = 0; k < data.size(); k++) {
- //****** coverage of AA **********//
- int index = -1;
- //find index in min where value is higher than d
- for (int m = 0; m < min1.size(); m++) {
- if (min1[m] > dist[k]) { index = m; break; }
- }
-
- // if you don't find one than all the mins are less than d
- if (index == -1) { index = min1.size(); }
-
- //****** coverage of AB **********//
- int index2 = -1;
- //find index in min where value is higher than d
- for (int m = 0; m < min2.size(); m++) {
- if (min2[m] > dist[k]) { index2 = m; break; }
- }
-
- // if you don't find one than all the mins are less than d
- if (index2 == -1) { index2 = min2.size(); }
-
- //coverage of ii
- float covII = 1.0 - ((min1.size()-index)/(float)min1.size());
-
- //coverage of ij
- float covIJ = 1.0 - ((min2.size()-index2)/(float)min2.size());
-
- //save value in data (Caa - Cab)^2 * distDiff
- data[k][i][j] = ((covII-covIJ) * (covII-covIJ)) * distDiff;
-
- //update distDiff
- if (k < data.size() - 1) {
- distDiff = dist[k+1] - dist[k];
- }
- }
-
- //put matrix back to original
- matrix->restore();
- }
- }
- count2++;
- }
- count += numGroups+1; //go from AA to BB to CC
- }
-
- }
- catch(exception& e) {
- cout << "Standard Error: " << e.what() << " has occurred in the Coverage class Function getValues. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
- exit(1);
- }
- catch(...) {
- cout << "An unknown error has occurred in the Coverage class function getValues. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
- exit(1);
- }
-
-}
-
+++ /dev/null
-#ifndef COVERAGE_H
-#define COVERAGE_H
-
-/*
- * coverage.h
- * Mothur
- *
- * Created by Sarah Westcott on 3/9/09.
- * Copyright 2009 Schloss Lab UMASS Amherst. All rights reserved.
- *
- */
-
-#include "mothur.h"
-#include "fullmatrix.h"
-#include "globaldata.hpp"
-
-using namespace std;
-
-/***********************************************************************/
-
-class Coverage {
-
- public:
- Coverage();
- ~Coverage(){};
- void getValues(FullMatrix*, vector< vector< vector<float> > >&, vector<float>, string); //matrix, container for results, vector of distances, mode - for random matrices
- void getValues(FullMatrix*, vector< vector< vector<float> > >&, vector<float>); //for user matrix
-
- private:
- GlobalData* globaldata;
- int numGroups, numUserGroups;
-
-};
-
-
-/***********************************************************************/
-
-
-
-#endif
\ No newline at end of file
--- /dev/null
+/*
+ * DLibshuff.cpp
+ * Mothur
+ *
+ * Created by Pat Schloss on 4/8/09.
+ * Copyright 2009 Patrick D. Schloss. All rights reserved.
+ *
+ */
+
+#include "dlibshuff.h"
+
+/***********************************************************************/
+
+DLibshuff::DLibshuff(FullMatrix* D, int it, float step, float co) : Libshuff(D, it, step, co){
+
+ numDXs = int(cutOff / stepSize);
+
+}
+
+/***********************************************************************/
+
+float DLibshuff::evaluatePair(int i, int j){
+ return dCalculate(i,j);
+}
+
+/***********************************************************************/
+
+vector<vector<double> > DLibshuff::evaluateAll(){
+ savedMins.resize(numGroups);
+ vector<vector<double> > dCXYValues(numGroups);
+
+ for(int i=0;i<numGroups;i++){
+ savedMins[i].resize(numGroups);
+ dCXYValues[i].resize(numGroups);
+ for(int j=0;j<numGroups;j++){
+ if(i!=j){ dCXYValues[i][j] = dCalculate(i,j); }
+ savedMins[i][i] = minX;
+ savedMins[i][j] = minXY;
+ }
+ }
+
+ return dCXYValues;
+}
+
+/***********************************************************************/
+
+double DLibshuff::dCalculate(int x, int y){
+
+ minX = getMinX(x);
+ minXY = getMinXY(x, y);
+
+ vector<int> nx = calcN(minX);
+ vector<int> nxy = calcN(minXY);
+
+ double sum = 0;
+
+ for(int i=0;i<numDXs;i++){
+ float h = (nx[i] - nxy[i]) / (float) groupSizes[x];
+ sum += h * h * stepSize;
+ }
+
+ return sum;
+}
+
+/***********************************************************************/
+
+vector<int> DLibshuff::calcN(vector<double> minVector){
+
+ vector<int> counts(numDXs,0);
+
+ int precision = int(1 / stepSize);
+
+ for(int i=0;i<minVector.size();i++){
+ int bin = int (precision * minVector[i]);
+ if(bin < numDXs){ counts[bin]++; }
+ }
+
+ for(int i=1;i<numDXs;i++){
+ counts[i] += counts[i-1];
+ }
+
+ return counts;
+}
+
+/***********************************************************************/
--- /dev/null
+#ifndef DLIBSHUFF
+#define DLIBSHUFF
+
+/*
+ * dlibshuff.h
+ * Mothur
+ *
+ * Created by Pat Schloss on 4/8/09.
+ * Copyright 2009 Patrick D. Schloss. All rights reserved.
+ *
+ */
+
+#include "fullmatrix.h"
+#include "libshuff.h"
+
+class DLibshuff : public Libshuff {
+
+public:
+ DLibshuff(FullMatrix*, int, float, float);
+ vector<vector<double> > evaluateAll();
+ float evaluatePair(int, int);
+
+private:
+ int numDXs;
+ double dCalculate(int, int);
+ vector<int> calcN(vector<double>);
+};
+
+#endif
#include "fullmatrix.h"
/**************************************************************************/
+
//This constructor reads a distance matrix file and stores the data in the matrix.
FullMatrix::FullMatrix(ifstream& filehandle) {
try{
group = groupmap->getGroup(name);
if(group == "not found") { cout << "Error: Sequence '" << name << "' was not found in the group file, please correct." << endl; exit(1); }
- index[0].groupname = group;
+ index.resize(numSeqs);
+ index[0].groupName = group;
index[0].seqName = name;
//determine if matrix is square or lower triangle
//if it is square read the distances for the first sequence
char d;
+ bool square;
while((d=filehandle.get()) != EOF){
//is d a number meaning its square
//sort sequences so they are gathered in groups for processing
sortGroups(0, numSeqs-1);
-
+
+ groups.push_back(index[0].groupName);
+ sizes.push_back(1);
+ int groupCount = 0;
+
+ for(int i=1;i<numSeqs;i++){
+ if(index[i].groupName == index[i-1].groupName){ sizes[groupCount]++; }
+ else{
+ sizes.push_back(1);
+ groups.push_back(index[i].groupName);
+ groupCount++;
+ }
+ }
+
}
catch(exception& e) {
cout << "Standard Error: " << e.what() << " has occurred in the FullMatrix class Function FullMatrix. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
filehandle >> name;
group = groupmap->getGroup(name);
- index[i].groupname = group;
+ index[i].groupName = group;
index[i].seqName = name;
if(group == "not found") { cout << "Error: Sequence '" << name << "' was not found in the group file, please correct." << endl; exit(1); }
filehandle >> name;
group = groupmap->getGroup(name);
- index[i].groupname = group;
+ index[i].groupName = group;
index[i].seqName = name;
if(group == "not found") { cout << "Error: Sequence '" << name << "' was not found in the group file, please correct." << endl; exit(1); }
}
/**************************************************************************/
+
void FullMatrix::sortGroups(int low, int high){
try{
/* compare value */
//what group does this row belong to
- string z = index[(low + high) / 2].groupname;
+ string z = index[(low + high) / 2].groupName;
/* partition */
do {
/* find member above ... */
- while(index[i].groupname < z) i++;
+ while(index[i].groupName < z) i++;
/* find element below ... */
- while(index[j].groupname > z) j--;
+ while(index[j].groupName > z) j--;
if(i <= j) {
/* swap rows*/
}
//swap map elements
- z = index[i].groupname;
- index[i].groupname = index[j].groupname;
- index[j].groupname = z;
+ z = index[i].groupName;
+ index[i].groupName = index[j].groupName;
+ index[j].groupName = z;
name = index[i].seqName;
index[i].seqName = index[j].seqName;
}
/**************************************************************************/
-int FullMatrix::getNumSeqs(){ return numSeqs; }
-/**************************************************************************/
-//print out matrix
-void FullMatrix::printMatrix(ostream& out) {
- try{
- for (int i = 0; i < numSeqs; i++) {
- out << "row " << i << " group = " << index[i].groupname << " name = " << index[i].seqName << endl;
- for (int j = 0; j < numSeqs; j++) {
- out << matrix[i][j] << " ";
- }
- out << endl;
- }
- }
- catch(exception& e) {
- cout << "Standard Error: " << e.what() << " has occurred in the FullMatrix class Function printMatrix. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
- exit(1);
- }
- catch(...) {
- cout << "An unknown error has occurred in the FullMatrix class function printMatrix. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
- exit(1);
- }
-
-}
-/**************************************************************************/
-void FullMatrix::setBounds(){
- try{
- numGroups = globaldata->gGroupmap->namesOfGroups.size();
-
- //sort globaldata->gGroupmap.namesOfGroups so that it will match the matrix
- sort(globaldata->gGroupmap->namesOfGroups.begin(), globaldata->gGroupmap->namesOfGroups.end());
-
- //one for each comparision
- //minsForRows.resize(numGroups*numGroups);
-
- /*************************************************/
- //find where in matrix each group starts and stops
- /*************************************************/
- bounds.resize(numGroups);
-
- bounds[0] = 0;
- bounds[numGroups] = numSeqs;
+float FullMatrix::get(int i, int j){ return matrix[i][j]; }
- //for each group find bounds of subgroup/comparison
- for (int i = 1; i < numGroups; i++) {
- getBounds(bounds[i], globaldata->gGroupmap->namesOfGroups[i-1]);
- }
-
- }
- catch(exception& e) {
- cout << "Standard Error: " << e.what() << " has occurred in the FullMatrix class Function getMinsForRowsVectors. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
- exit(1);
- }
- catch(...) {
- cout << "An unknown error has occurred in the FullMatrix class function getMinsForRowsVectors. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
- exit(1);
- }
+/**************************************************************************/
-}
-/**************************************************************************/
-vector<float> FullMatrix::getMins(int x) {
- try{
- //clear out old data
- minsForRows.clear();
-
- /************************************************************/
- //fill the minsForRows vector for the box the user wants
- /************************************************************/
- int count = 0;
- int lowBoundx = bounds[0]; //where first group starts
- int lowBoundy = bounds[0];
- int highBoundx = bounds[1]; //where second group starts
- int highBoundy = bounds[1];
-
- int countx = 1; //index in bound
- int county = 1; //index in bound
-
- //find the bounds for the box the user wants
- for (int i = 0; i < (numGroups * numGroups); i++) {
-
- //are you at the box?
- if (count == x) { break; }
- else { count++; }
-
- //move to next box
- if (county < numGroups) {
- county++;
- highBoundy = bounds[county];
- lowBoundy = bounds[county-1];
- }else{ //you are moving to a new row of "boxes"
- county = 1;
- countx++;
- highBoundx = bounds[countx];
- lowBoundx = bounds[countx-1];
- highBoundy = bounds[county];
- lowBoundy = bounds[county-1];
- }
- }
-
- //each row in the box
- for (int x = lowBoundx; x < highBoundx; x++) {
- float min4Row = 100000.0;
- //each entry in that row
- for (int y = lowBoundy; y < highBoundy; y++) {
- //if you are not on the diagonal and you are less than previous minimum
- if ((x != y) && (matrix[x][y] < min4Row)) {
- min4Row = matrix[x][y];
- }
- }
- //save minimum value for that row in minsForRows vector of vectors
- minsForRows.push_back(min4Row);
- }
-
- return minsForRows;
- }
- catch(exception& e) {
- cout << "Standard Error: " << e.what() << " has occurred in the FullMatrix class Function getMins. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
- exit(1);
- }
- catch(...) {
- cout << "An unknown error has occurred in the FullMatrix class function getMins. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
- exit(1);
- }
-}
-/**************************************************************************/
-void FullMatrix::getBounds(int& higher, string group) {
- try{
- bool gotLower = false;
-
- //for each group find bounds of subgroup/comparison
- for (it = index.begin(); it != index.end(); it++) {
- if (it->second.groupname == group) {
- gotLower = true;
- }else if ((gotLower == true) && (it->second.groupname != group)) { higher = it->first; break; }
- }
-
- }
- catch(exception& e) {
- cout << "Standard Error: " << e.what() << " has occurred in the FullMatrix class Function getBounds. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
- exit(1);
- }
- catch(...) {
- cout << "An unknown error has occurred in the FullMatrix class function getBounds. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
- exit(1);
- }
+vector<string> FullMatrix::getGroups(){ return groups; }
-}
+/**************************************************************************/
-/**************************************************************************/
-//print out matrix
-void FullMatrix::printMinsForRows(ostream& out) {
- try{
- for (int j = 0; j < minsForRows.size(); j++) {
- out << minsForRows[j] << " ";
- }
- out << endl;
+vector<int> FullMatrix::getSizes(){ return sizes; }
- }
- catch(exception& e) {
- cout << "Standard Error: " << e.what() << " has occurred in the FullMatrix class Function printMatrix. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
- exit(1);
- }
- catch(...) {
- cout << "An unknown error has occurred in the FullMatrix class function printMatrix. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
- exit(1);
- }
+/**************************************************************************/
-}
+int FullMatrix::getNumGroups(){ return groups.size(); }
-/**************************************************************************/
-//shuffles the sequences in the 2 groups passed in.
-void FullMatrix::shuffle(string groupA, string groupB){
- try{
- vector<int> rows2Swap;
- vector<int> shuffled;
- float y = 0;
- string name = "";
-
-
- /********************************/
- //save rows you want to randomize
- /********************************/
- //go through the matrix map to find the rows from groups you want to randomize
- for (it = index.begin(); it != index.end(); it++) {
- //is this row from group A or B?
- if ((it->second.groupname == groupA) || (it->second.groupname == groupB)) {
- rows2Swap.push_back(it->first);
- shuffled.push_back(it->first);
- }
- }
-
- //randomize rows to shuffle in shuffled
- random_shuffle(shuffled.begin(), shuffled.end());
-
- /***************************************/
- //swap rows and columns to randomize box
- /***************************************/
- for (int i = 0; i < shuffled.size(); i++) {
+/**************************************************************************/
- //record the swaps you are making so you can undo them in restore function
- restoreIndex[i].a = shuffled[i];
- restoreIndex[i].b = rows2Swap[i];
-
- /* swap rows*/
- for (int h = 0; h < numSeqs; h++) {
- y = matrix[shuffled[i]][h];
- matrix[shuffled[i]][h] = matrix[rows2Swap[i]][h];
- matrix[rows2Swap[i]][h] = y;
- }
-
- /* swap columns */
- for (int b = 0; b < numSeqs; b++) {
- y = matrix[b][shuffled[i]];
- matrix[b][shuffled[i]] = matrix[b][rows2Swap[i]];
- matrix[b][rows2Swap[i]] = y;
- }
-
- //swap map elements
- name = index[shuffled[i]].seqName;
- index[shuffled[i]].seqName = index[rows2Swap[i]].seqName;
- index[rows2Swap[i]].seqName = name;
+int FullMatrix::getNumSeqs(){ return numSeqs; }
- }
- }
- catch(exception& e) {
- cout << "Standard Error: " << e.what() << " has occurred in the FullMatrix class Function shuffle. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
- exit(1);
- }
- catch(...) {
- cout << "An unknown error has occurred in the FullMatrix class function shuffle. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
- exit(1);
- }
-}
/**************************************************************************/
-//unshuffles the matrix.
-void FullMatrix::restore(){
- try{
- float y = 0;
- string name = "";
-
- //reverse iterate through swaps and undo them to restore original matrix and index map.
- for(it2 = restoreIndex.rbegin(); it2 != restoreIndex.rend(); it2++) {
- /* swap rows */
-
- for (int h = 0; h < numSeqs; h++) {
- y = matrix[it2->second.a][h];
- matrix[it2->second.a][h] = matrix[it2->second.b][h];
- matrix[it2->second.b][h] = y;
- }
-
- /* swap columns */
- for (int b = 0; b < numSeqs; b++) {
- y = matrix[b][it2->second.a];
- matrix[b][it2->second.a] = matrix[b][it2->second.b];
- matrix[b][it2->second.b] = y;
- }
-
-
- //swap map elements
- name = index[it2->second.a].seqName;
- index[it2->second.a].seqName = index[it2->second.b].seqName;
- index[it2->second.b].seqName = name;
- }
-
- //clear restore for next shuffle
- restoreIndex.clear();
- }
- catch(exception& e) {
- cout << "Standard Error: " << e.what() << " has occurred in the FullMatrix class Function restore. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
- exit(1);
- }
- catch(...) {
- cout << "An unknown error has occurred in the FullMatrix class function restore. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
- exit(1);
- }
-}
-/**************************************************************************/
-void FullMatrix::getDist(vector<float>& distances) {
+void FullMatrix::printMatrix(ostream& out) {
try{
- map<float, float> dist; //holds the distances for the integral form
- map<float, float>::iterator it;
-
- /************************************************************/
- //fill the minsForRows vectors for each group the user wants
- /************************************************************/
- int lowBoundx = bounds[0]; //where first group starts
- int lowBoundy = bounds[0];
- int highBoundx = bounds[1]; //where second group starts
- int highBoundy = bounds[1];
-
- int countx = 1; //index in bound
- int county = 1; //index in bound
-
- //go through each "box" in the matrix
- for (int i = 0; i < (numGroups * numGroups); i++) {
- //each row in the box
- for (int x = lowBoundx; x < highBoundx; x++) {
- float min4Row = 100000.0;
- //each entry in that row
- for (int y = lowBoundy; y < highBoundy; y++) {
- //if you are not on the diagonal and you are less than previous minimum
- if ((x != y) && (matrix[x][y] < min4Row)){
- min4Row = matrix[x][y];
- }
- }
- //save minimum value
- dist[min4Row] = min4Row;
- }
-
- //****** reset bounds to process next "box" ********
- //if you still have more "boxes" in that row
- if (county < numGroups) {
- county++;
- highBoundy = bounds[county];
- lowBoundy = bounds[county-1];
- }else{ //you are moving to a new row of "boxes"
- county = 1;
- countx++;
- highBoundx = bounds[countx];
- lowBoundx = bounds[countx-1];
- highBoundy = bounds[county];
- lowBoundy = bounds[county-1];
+ for (int i = 0; i < numSeqs; i++) {
+ out << "row " << i << " group = " << index[i].groupName << " name = " << index[i].seqName << endl;
+ for (int j = 0; j < numSeqs; j++) {
+ out << matrix[i][j] << " ";
}
+ out << endl;
}
-
- //store distances in users vector
- for (it = dist.begin(); it != dist.end(); it++) {
- distances.push_back(it->first);
- }
-
}
catch(exception& e) {
- cout << "Standard Error: " << e.what() << " has occurred in the FullMatrix class Function restore. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
+ cout << "Standard Error: " << e.what() << " has occurred in the FullMatrix class Function printMatrix. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
exit(1);
}
catch(...) {
- cout << "An unknown error has occurred in the FullMatrix class function restore. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
+ cout << "An unknown error has occurred in the FullMatrix class function printMatrix. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
exit(1);
}
+
}
+/**************************************************************************/
+
using namespace std;
struct Names {
- string groupname;
string seqName;
+ string groupName;
};
-struct Swap {
- int a;
- int b;
-};
-
-
class FullMatrix {
- public:
- FullMatrix(){};
- FullMatrix(ifstream&);
- ~FullMatrix(){};
+public:
+ FullMatrix(){};
+ FullMatrix(ifstream&);
+ ~FullMatrix(){};
- int getNumSeqs();
- void printMatrix(ostream&);
- void setBounds(); //requires globaldata->Groups to be filled
- vector<float> getMins(int); //returns vector of mins for "box" requested ie. groups A, B, 0 = AA, 1 = AB, 2 = BA, 3 = BB;
- void getDist(vector<float>&); //fills a vector with the valid distances for the integral form.
- void shuffle(string, string); //shuffles the sequences in the groups passed in.
- void restore(); //unshuffles the matrix.
+ int getNumSeqs();
+ vector<int> getSizes();
+ vector<string> getGroups();
+ int getNumGroups();
+ void printMatrix(ostream&);
+ float get(int, int);
- private:
- void sortGroups(int, int); //this function sorts the sequences within the matrix.
- void getBounds(int&, string);
- void readSquareMatrix(ifstream&);
- void readLTMatrix(ifstream&);
- void printMinsForRows(ostream&);
-
- map<int, Names> index; // row in vector, sequence group. need to know this so when we sort it can be updated.
- map<int, Swap> restoreIndex; //a map of the swaps made so you can undo them in restore.
- map<int, Names>::iterator it;
- map<int, Swap>::reverse_iterator it2;
-
- vector< vector<float> > matrix; //a 2D distance matrix of all the sequences and their distances to eachother.
- vector<float> minsForRows; //vector< minimum distance for that subrow> - one for each comparison.
- vector<int> bounds; //bounds[1] = starting row in matrix from group B, bounds[2] = starting row in matrix from group C, bounds[3] = no need to find upper bound of C because its numSeqs.
-
- GroupMap* groupmap; //maps sequences to groups they belong to.
- GlobalData* globaldata;
- int numSeqs, numGroups, numUserGroups;
- bool square;
-
+private:
+ vector< vector<float> > matrix; //a 2D distance matrix of all the sequences and their distances to eachother.
+ void readSquareMatrix(ifstream&);
+ void readLTMatrix(ifstream&);
+ vector<Names> index; // row in vector, sequence group. need to know this so when we sort it can be updated.
+ vector<int> sizes;
+ vector<string> groups;
+ void sortGroups(int, int); //this function sorts the sequences within the matrix.
+
+
+ GroupMap* groupmap; //maps sequences to groups they belong to.
+ int numSeqs;
+ int numGroups;
+ int numUserGroups;
+ GlobalData* globaldata;
};
#endif
\ No newline at end of file
--- /dev/null
+/*
+ * libshuffform.cpp
+ * Mothur
+ *
+ * Created by Pat Schloss on 4/8/09.
+ * Copyright 2009 Patrick D. Schloss. All rights reserved.
+ *
+ */
+
+#include "libshuff.h"
+
+/***********************************************************************/
+
+void swap(int& i,int& j){ int t = i; i = j; j = t; }
+
+/***********************************************************************/
+
+Libshuff::Libshuff(FullMatrix* D, int it, float step, float co) : matrix(D), iters(it), stepSize(step), cutOff(co){
+ try{
+ groupNames = matrix->getGroups();
+ groupSizes = matrix->getSizes();
+ numGroups = matrix->getNumGroups();
+
+ initializeGroups(matrix);
+ }
+ catch(exception& e) {
+ cout << "Standard Error: " << e.what() << " has occurred in the Libshuff class Function Libshuff. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
+ exit(1);
+ }
+ catch(...) {
+ cout << "An unknown error has occurred in the Libshuff class function Libshuff. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
+ exit(1);
+ }
+
+}
+
+/***********************************************************************/
+
+void Libshuff::initializeGroups(FullMatrix* matrix){
+ try{
+ groups.resize(numGroups);
+ savedGroups.resize(numGroups);
+
+ savedGroups.resize(numGroups);
+ for(int i=0;i<numGroups;i++) {
+ groups[i].resize(groupSizes[i]);
+ savedGroups[i].resize(groupSizes[i]);
+ }
+ int index=0;
+ for(int i=0;i<numGroups;i++){
+ for(int j=0;j<groupSizes[i];j++){
+ savedGroups[i][j] = groups[i][j] = index++;
+ }
+ }
+ }
+ catch(exception& e) {
+ cout << "Standard Error: " << e.what() << " has occurred in the Libshuff class Function initializeGroups. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
+ exit(1);
+ }
+ catch(...) {
+ cout << "An unknown error has occurred in the Libshuff class function initializeGroups. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
+ exit(1);
+ }
+
+}
+
+/***********************************************************************/
+
+vector<vector<vector<double> > > Libshuff::getSavedMins(){
+ return savedMins;
+}
+
+/***********************************************************************/
+
+vector<double> Libshuff::getMinX(int x){
+ try{
+ vector<double> minX(groupSizes[x], 0);
+ for(int i=0;i<groupSizes[x];i++){
+ minX[i] = (groupSizes[x] > 1 ? (i==0 ? matrix->get(groups[x][0], groups[x][1]) : matrix->get(groups[x][i], groups[x][0])) : 0.0);
+ for(int j=0;j<groupSizes[x];j++){
+ if(i != j) {
+ double dx = matrix->get(groups[x][i], groups[x][j]);
+ if(dx < minX[i]){ minX[i] = dx; }
+ }
+ }
+ }
+ return minX;
+ }
+ catch(exception& e) {
+ cout << "Standard Error: " << e.what() << " has occurred in the Libshuff class Function getMinX. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
+ exit(1);
+ }
+ catch(...) {
+ cout << "An unknown error has occurred in the Libshuff class function getMinX. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
+ exit(1);
+ }
+
+}
+
+/***********************************************************************/
+
+vector<double> Libshuff::getMinXY(int x, int y){
+ try{
+ vector<double> minXY(groupSizes[x], 0);
+
+ for(int i=0;i<groupSizes[x];i++){
+ minXY[i] = matrix->get(groups[x][i], groups[y][0]);
+ for(int j=0;j<groupSizes[y];j++){
+ double dxy = matrix->get(groups[x][i], groups[y][j]);
+ if(dxy<minXY[i]){ minXY[i] = dxy; }
+ }
+ }
+ return minXY;
+ }
+ catch(exception& e) {
+ cout << "Standard Error: " << e.what() << " has occurred in the Libshuff class Function getMinXY. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
+ exit(1);
+ }
+ catch(...) {
+ cout << "An unknown error has occurred in the Libshuff class function getMinXY. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
+ exit(1);
+ }
+}
+
+/***********************************************************************/
+
+void Libshuff::randomizeGroups(int x, int y){
+ try{
+ int nv = groupSizes[x]+groupSizes[y];
+ vector<int> v(nv);
+
+ int index=0;
+ for(int k=0;k<groupSizes[x];k++) { v[index++] = groups[x][k]; }
+ for(int k=0;k<groupSizes[y];k++) { v[index++] = groups[y][k]; }
+
+ for(int k=nv-1;k>0;k--){
+ int z = (int)(rand() % k);
+ swap(v[z],v[k]);
+ }
+
+ index=0;
+ for(int k=0;k<groupSizes[x];k++) { groups[x][k]=v[index++]; }
+ for(int k=0;k<groupSizes[y];k++) { groups[y][k]=v[index++]; }
+ }
+ catch(exception& e) {
+ cout << "Standard Error: " << e.what() << " has occurred in the Libshuff class Function randomizeGroups. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
+ exit(1);
+ }
+ catch(...) {
+ cout << "An unknown error has occurred in the Libshuff class function randomizeGroups. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
+ exit(1);
+ }
+}
+
+/***********************************************************************/
+
+void Libshuff::resetGroup(int x){
+
+ for(int k=0;k<groupSizes[x];k++) { groups[x][k] = savedGroups[x][k]; }
+
+}
+
+/***********************************************************************/
--- /dev/null
+#ifndef LIBSHUFF
+#define LIBSHUFF
+
+/*
+ * libshuff.h
+ * Mothur
+ *
+ * Created by Pat Schloss on 4/8/09.
+ * Copyright 2009 Patrick D. Schloss. All rights reserved.
+ *
+ */
+
+#include "fullmatrix.h"
+
+class Libshuff {
+
+public:
+ Libshuff(FullMatrix*, int, float, float);
+ virtual vector<vector<double> > evaluateAll() = 0;
+ virtual float evaluatePair(int,int) = 0;
+ void randomizeGroups(int, int);
+ void resetGroup(int);
+ vector<vector<vector<double> > > getSavedMins();
+
+protected:
+ void initializeGroups(FullMatrix*);
+ vector<double> getMinX(int);
+ vector<double> getMinXY(int, int);
+
+ vector<vector<vector<double> > > savedMins;
+
+
+ FullMatrix* matrix;
+ vector<int> groupSizes;
+ vector<string> groupNames;
+ vector<vector<int> > groups;
+ vector<vector<int> > savedGroups;
+ vector<double> minX;
+ vector<double> minXY;
+ float cutOff;
+ int iters;
+ float stepSize;
+
+ int numGroups;
+};
+
+#endif
#include "libshuffcommand.h"
+#include "libshuff.h"
+#include "slibshuff.h"
+#include "dlibshuff.h"
//**********************************************************************************************************************
-
LibShuffCommand::LibShuffCommand(){
try {
- globaldata = GlobalData::getInstance();
- convert(globaldata->getCutOff(), cutOff);
- convert(globaldata->getIters(), iters);
- convert(globaldata->getStep(), step);
- form = globaldata->getForm();
- matrix = globaldata->gMatrix;
- coverageFile = getRootName(globaldata->getPhylipFile()) + "coverage";
- summaryFile = getRootName(globaldata->getPhylipFile()) + "slsummary";
- openOutputFile(coverageFile, out);
- openOutputFile(summaryFile, outSum);
+ srand( (unsigned)time( NULL ) );
- //set the groups to be analyzed
- setGroups();
+ globaldata = GlobalData::getInstance();
+ convert(globaldata->getCutOff(), cutOff); //get the cutoff
+ convert(globaldata->getIters(), iters); //get the number of iterations
+ convert(globaldata->getStep(), step); //get the step size for the discrete command
+ matrix = globaldata->gMatrix; //get the distance matrix
+ setGroups(); //set the groups to be analyzed
- //file headers for coverage file
- out << "D" << '\t';
- for (int i = 0; i < groupComb.size(); i++) {
- out << "C" + groupComb[i] << '\t';
+ if(globaldata->getForm() == "discrete"){
+ form = new DLibshuff(matrix, iters, step, cutOff);
}
-
- for (int i = 0; i < numGroups; i++) {
- for (int j = 0; j < numGroups; j++) {
- //don't output AA to AA
- if (i != j) {
- out << "Delta" + globaldata->Groups[i] + "-" + globaldata->Groups[j] << '\t';
- }
- }
+ else{
+ form = new SLibshuff(matrix, iters, cutOff);
}
- out << endl;
-
- numComp = numGroups*numGroups;
-
- coverage = new Coverage();
-
}
catch(exception& e) {
cout << "Standard Error: " << e.what() << " has occurred in the LibShuffCommand class Function LibShuffCommand. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
//**********************************************************************************************************************
-LibShuffCommand::~LibShuffCommand(){
- delete coverage;
-}
-
-//**********************************************************************************************************************
-
int LibShuffCommand::execute(){
try {
- //deltaValues[0] = scores for the difference between AA and AB.
- //cValues[0][0][0] = AA at distance 0.0, cValues[0][0][1] = AB at distance 0.0, cValues[0][0][2] = AC at distance 0.0, cValues[0][1][0] = BA at distance 0.0, cValues[0][1][1] = BB...
- Progress* reading;
- reading = new Progress("Comparing to random:", iters);
-
- sumDelta.resize(numComp-numGroups, 0.0);
-
- matrix->setBounds();
-
- //load distances
- if (form != "discrete") { matrix->getDist(dist); }
- else {
- float f = 0.0;
- while (f <= cutOff) {
- dist.push_back(f);
- f += step;
- }
- }
-
- /*****************************/
- //get values for users matrix
- /*****************************/
-
- //clear out old Values
- deltaValues.clear();
- deltaValues.resize(dist.size());
-
- coverage->getValues(matrix, cValues, dist);
-
- float distDiff = 0;
-
- //loop through each distance and load sumdelta
- for (int p = 0; p < cValues.size(); p++) {
- //find delta values
- int count = 0;
- for (int i = 0; i < numGroups; i++) {
- for (int j = 0; j < numGroups; j++) {
- //don't save AA to AA
- if (i != j) {
- //(Caa - Cab)^2 * distDiff
- deltaValues[p].push_back(((cValues[p][i][i]-cValues[p][i][j]) * (cValues[p][i][i]-cValues[p][i][j])) * distDiff); //* distDiff
- sumDelta[count] += deltaValues[p][count];
- count++;
- }
- }
- }
- if (p < cValues.size() - 1) {
- distDiff = dist[p+1] - dist[p];
- }
- }
-
- printCoverageFile();
-
- /*******************************************************************************/
- //create and score random matrixes finding the sumDelta values for summary file
- /******************************************************************************/
- //initialize rsumDelta
- rsumDelta.resize(numComp-numGroups);
- for (int l = 0; l < rsumDelta.size(); l++) {
- for (int w = 0; w < iters; w++) {
- rsumDelta[l].push_back(0.0);
- }
- }
+ savedDXYValues = form->evaluateAll();
+ savedMinValues = form->getSavedMins();
+ pValueCounts.assign(numGroups, 0);
+ for(int i=0;i<numGroups;i++){
+ pValueCounts[i].assign(numGroups, 0);
+ }
- for (int m = 0; m < iters; m++) {
- //generate random matrix in getValues
- //values for random matrix
+ Progress* reading = new Progress();
- coverage->getValues(matrix, cValues, dist, "random");
-
- //loop through each distance and load rsumdelta
- for (int p = 0; p < cValues.size(); p++) {
- //find delta values
- int count = 0;
- for (int i = 0; i < numGroups; i++) {
- for (int j = 0; j < numGroups; j++) {
- //don't save AA to AA
- if (i != j) {
- //rsumDelta[3][500] = the sum of the delta scores for BB-BC for random matrix # 500.
- rsumDelta[count][m] += cValues[p][i][j]; // where cValues[p][0][1] = delta value at distance p of AA-AB, cValues[p][1][2] = delta value at distance p of BB-BC.
- count++;
- }
- }
+ for(int i=0;i<numGroups-1;i++) {
+ for(int j=i+1;j<numGroups;j++) {
+ reading->newLine(groupNames[i]+'-'+groupNames[j], iters);
+ for(int p=0;p<iters;p++) {
+ form->randomizeGroups(i,j);
+ if(form->evaluatePair(i,j) >= savedDXYValues[i][j]) { pValueCounts[i][j]++; }
+ if(form->evaluatePair(j,i) >= savedDXYValues[j][i]) { pValueCounts[j][i]++; }
+ reading->update(p);
}
+ form->resetGroup(i);
+ form->resetGroup(j);
}
-
- //clear out old Values
- reading->update(m);
- cValues.clear();
-
}
-
reading->finish();
delete reading;
-
- /**********************************************************/
- //find the signifigance of the user matrix' sumdelta values
- /**********************************************************/
-
- for (int t = 0; t < rsumDelta.size(); t++) {
- //sort rsumDelta t
- sort(rsumDelta[t].begin(), rsumDelta[t].end());
-
- //the index of the score higher than yours is returned
- //so if you have 1000 random matrices the index returned is 100
- //then there are 900 matrices with a score greater then you.
- //giving you a signifigance of 0.900
- int index = findIndex(sumDelta[t], t);
-
- //the signifigance is the number of trees with the users score or higher
- sumDeltaSig.push_back((iters-index)/(float)iters);
- }
-
+ cout << endl;
printSummaryFile();
+ printCoverageFile();
//clear out users groups
globaldata->Groups.clear();
+ delete form;
return 0;
}
exit(1);
}
}
+
//**********************************************************************************************************************
+
void LibShuffCommand::printCoverageFile() {
try {
- //format output
- out.setf(ios::fixed, ios::floatfield); out.setf(ios::showpoint);
+
+ ofstream outCov;
+ summaryFile = getRootName(globaldata->getPhylipFile()) + "libshuff.coverage";
+ openOutputFile(summaryFile, outCov);
+ outCov.setf(ios::fixed, ios::floatfield); outCov.setf(ios::showpoint);
+ cout.setf(ios::fixed, ios::floatfield); cout.setf(ios::showpoint);
+
+ map<double,vector<int> > allDistances;
+ map<double,vector<int> >::iterator it;
+
+ vector<vector<int> > indices(numGroups);
+ int numIndices = numGroups * numGroups;
- //loop through each distance
- for (int p = 0; p < cValues.size(); p++) {
- out << setprecision(6) << dist[p] << '\t';
- //print out coverage values
- for (int i = 0; i < numGroups; i++) {
- for (int j = 0; j < numGroups; j++) {
- out << cValues[p][i][j] << '\t';
+ int index = 0;
+ for(int i=0;i<numGroups;i++){
+ indices[i].assign(numGroups,0);
+ for(int j=0;j<numGroups;j++){
+ indices[i][j] = index++;
+ for(int k=0;k<savedMinValues[i][j].size();k++){
+ if(allDistances[savedMinValues[i][j][k]].size() != 0){
+ allDistances[savedMinValues[i][j][k]][indices[i][j]]++;
+ }
+ else{
+ allDistances[savedMinValues[i][j][k]].assign(numIndices, 0);
+ allDistances[savedMinValues[i][j][k]][indices[i][j]] = 1;
+ }
}
}
-
- for (int h = 0; h < deltaValues[p].size(); h++) {
- out << deltaValues[p][h] << '\t';
+ }
+ it=allDistances.begin();
+
+ cout << setprecision(8);
+
+ vector<int> prevRow = it->second;
+ it++;
+
+ for(it;it!=allDistances.end();it++){
+ for(int i=0;i<it->second.size();i++){
+ it->second[i] += prevRow[i];
+ }
+ prevRow = it->second;
+ }
+
+ vector<int> lastRow = allDistances.rbegin()->second;
+ outCov << setprecision(8);
+
+ outCov << "dist";
+ for (int i = 0; i < numGroups; i++){
+ outCov << '\t' << groupNames[i];
+ }
+ for (int i=0;i<numGroups;i++){
+ for(int j=i+1;j<numGroups;j++){
+ outCov << '\t' << groupNames[i] << '-' << groupNames[j] << '\t';
+ outCov << groupNames[j] << '-' << groupNames[i];
}
-
- out << endl;
+ }
+ outCov << endl;
+
+ for(it=allDistances.begin();it!=allDistances.end();it++){
+ outCov << it->first << '\t';
+ for(int i=0;i<numGroups;i++){
+ outCov << it->second[indices[i][i]]/(float)lastRow[indices[i][i]] << '\t';
+ }
+ for(int i=0;i<numGroups;i++){
+ for(int j=i+1;j<numGroups;j++){
+ outCov << it->second[indices[i][j]]/(float)lastRow[indices[i][j]] << '\t';
+ outCov << it->second[indices[j][i]]/(float)lastRow[indices[j][i]] << '\t';
+ }
+ }
+ outCov << endl;
}
}
exit(1);
}
}
+
//**********************************************************************************************************************
+
void LibShuffCommand::printSummaryFile() {
try {
- //format output
+
+ ofstream outSum;
+ summaryFile = getRootName(globaldata->getPhylipFile()) + "libshuff.summary";
+ openOutputFile(summaryFile, outSum);
+
outSum.setf(ios::fixed, ios::floatfield); outSum.setf(ios::showpoint);
+ cout.setf(ios::fixed, ios::floatfield); cout.setf(ios::showpoint);
- for (int i = 0; i < numGroups; i++) {
- for (int j = 0; j < numGroups; j++) {
- //don't output AA to AA
- if (i != j) {
- outSum << "Delta " + globaldata->Groups[i] + "-" + globaldata->Groups[j] << '\t'<< "DeltaSig " + globaldata->Groups[i] + "-" + globaldata->Groups[j] << '\t';
- cout << "Delta " + globaldata->Groups[i] + "-" + globaldata->Groups[j] << '\t'<< "DeltaSig " + globaldata->Groups[i] + "-" + globaldata->Groups[j] << '\t';
+ cout << setw(20) << left << "Comparison" << '\t' << setprecision(8) << "dCXYScore" << '\t' << "Significance" << endl;
+ outSum << setw(20) << left << "Comparison" << '\t' << setprecision(8) << "dCXYScore" << '\t' << "Significance" << endl;
+
+ int precision = (int)log10(iters);
+ for(int i=0;i<numGroups;i++){
+ for(int j=i+1;j<numGroups;j++){
+ if(pValueCounts[i][j]){
+ cout << setw(20) << left << groupNames[i]+'-'+groupNames[j] << '\t' << setprecision(8) << savedDXYValues[i][j] << '\t' << setprecision(precision) << pValueCounts[i][j]/(float)iters << endl;
+ outSum << setw(20) << left << groupNames[i]+'-'+groupNames[j] << '\t' << setprecision(8) << savedDXYValues[i][j] << '\t' << setprecision(precision) << pValueCounts[i][j]/(float)iters << endl;
+ }
+ else{
+ cout << setw(20) << left << groupNames[i]+'-'+groupNames[j] << '\t' << setprecision(8) << savedDXYValues[i][j] << '\t' << '<' <<setprecision(precision) << 1/(float)iters << endl;
+ outSum << setw(20) << left << groupNames[i]+'-'+groupNames[j] << '\t' << setprecision(8) << savedDXYValues[i][j] << '\t' << '<' <<setprecision(precision) << 1/(float)iters << endl;
+ }
+ if(pValueCounts[j][i]){
+ cout << setw(20) << left << groupNames[j]+'-'+groupNames[i] << '\t' << setprecision(8) << savedDXYValues[j][i] << '\t' << setprecision (precision) << pValueCounts[j][i]/(float)iters << endl;
+ outSum << setw(20) << left << groupNames[j]+'-'+groupNames[i] << '\t' << setprecision(8) << savedDXYValues[j][i] << '\t' << setprecision (precision) << pValueCounts[j][i]/(float)iters << endl;
+ }
+ else{
+ cout << setw(20) << left << groupNames[j]+'-'+groupNames[i] << '\t' << setprecision(8) << savedDXYValues[j][i] << '\t' << '<' <<setprecision (precision) << 1/(float)iters << endl;
+ outSum << setw(20) << left << groupNames[j]+'-'+groupNames[i] << '\t' << setprecision(8) << savedDXYValues[j][i] << '\t' << '<' <<setprecision (precision) << 1/(float)iters << endl;
}
}
}
- outSum << endl;
- cout << endl;
- //print out delta values
- for (int i = 0; i < sumDelta.size(); i++) {
- if (sumDeltaSig[i] > (1/(float)iters)) {
- outSum << setprecision(6) << sumDelta[i] << '\t' << setprecision(globaldata->getIters().length()) << sumDeltaSig[i] << '\t';
- cout << setprecision(6) << sumDelta[i] << '\t' << setprecision(globaldata->getIters().length()) << sumDeltaSig[i] << '\t';
- }else {
- outSum << setprecision(6) << sumDelta[i] << '\t' << setprecision(globaldata->getIters().length()) << "<" << (1/float(iters)) << '\t';
- cout << setprecision(6) << sumDelta[i] << '\t' << setprecision(globaldata->getIters().length()) << "<" << (1/float(iters)) << '\t';
- }
- }
- outSum << endl;
- cout << endl;
-
}
catch(exception& e) {
cout << "Standard Error: " << e.what() << " has occurred in the LibShuffCommand class Function printSummaryFile. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
}
//**********************************************************************************************************************
+
void LibShuffCommand::setGroups() {
try {
//if the user has not entered specific groups to analyze then do them all
for (int i=0; i < numGroups; i++) {
globaldata->Groups.push_back(globaldata->gGroupmap->namesOfGroups[i]);
}
- }else {
+ } else {
if (globaldata->getGroups() != "all") {
//check that groups are valid
for (int i = 0; i < globaldata->Groups.size(); i++) {
globaldata->Groups.push_back(globaldata->gGroupmap->namesOfGroups[i]);
}
cout << "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." << endl;
- }else { numGroups = globaldata->Groups.size(); }
- }else { //users wants all groups
+ } else { numGroups = globaldata->Groups.size(); }
+ } else { //users wants all groups
numGroups = globaldata->gGroupmap->getNumGroups();
globaldata->Groups.clear();
for (int i=0; i < numGroups; i++) {
}
}
}
-
+
//sort so labels match
sort(globaldata->Groups.begin(), globaldata->Groups.end());
//sort
sort(globaldata->gGroupmap->namesOfGroups.begin(), globaldata->gGroupmap->namesOfGroups.end());
-
+ groupNames = globaldata->Groups;
+
// number of comparisons i.e. with groups A,B,C = AA, AB, AC, BA, BB, BC...;
- for (int i=0; i<numGroups; i++) {
- for (int l = 0; l < numGroups; l++) {
- //set group comparison labels
- groupComb.push_back(globaldata->Groups[i] + "-" + globaldata->Groups[l]);
- }
- }
+// for (int i=0; i<numGroups; i++) {
+// for (int l = 0; l < numGroups; l++) {
+// //set group comparison labels
+// groupComb.push_back(globaldata->Groups[i] + "-" + globaldata->Groups[l]);
+// }
+// }
}
catch(exception& e) {
cout << "Standard Error: " << e.what() << " has occurred in the LibShuffCommand class Function setGroups. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
exit(1);
}
}
-/***********************************************************/
-int LibShuffCommand::findIndex(float score, int index) {
- try{
- for (int i = 0; i < rsumDelta[index].size(); i++) {
- if (rsumDelta[index][i] >= score) { return i; }
- }
- return rsumDelta[index].size();
- }
- catch(exception& e) {
- cout << "Standard Error: " << e.what() << " has occurred in the LibShuffCommand class Function findIndex. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
- exit(1);
- }
- catch(...) {
- cout << "An unknown error has occurred in the LibShuffCommand class function findIndex. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
- exit(1);
- }
-}
/***********************************************************/
-
*/
#include "command.hpp"
-#include "coverage.h"
#include "fullmatrix.h"
+#include "libshuff.h"
using namespace std;
public:
LibShuffCommand();
- ~LibShuffCommand();
+ ~LibShuffCommand(){};
int execute();
private:
- vector< vector< vector<float> > > cValues; // vector<vector of coverage scores, one for each comparison.> -one for each distance level.
- vector< vector<float> > deltaValues; // vector< vector of delta scores, one for each comparison.> -one at each distance
- vector<float> sumDelta; //sum of delta scores, one for each comparison.
- vector<float> sumDeltaSig; //number of random matrixes with that delta value or ??
- vector< vector<float> > rsumDelta; //vector< vector<sumdelta scores for a given comparison> >
- vector<string> groupComb;
- vector<float> dist;
-
+ vector<string> groupNames;
void setGroups();
- int findIndex(float, int);
void printCoverageFile();
void printSummaryFile();
GlobalData* globaldata;
- Coverage* coverage;
FullMatrix* matrix;
+ Libshuff* form;
float cutOff, step;
int numGroups, numComp, iters;
- string coverageFile, summaryFile, form;
- ofstream out, outSum;
-
-
+ string coverageFile, summaryFile;
+ vector<vector<int> > pValueCounts;
+ vector<vector<double> > savedDXYValues;
+ vector<vector<vector<double> > > savedMinValues;
};
#endif
const char marker = '|';
+/***********************************************************************/
+
+Progress::Progress(){
+ try {
+ cout << "********************#****#****#****#****#****#****#****#****#****#****#";
+
+ nTicks = 0;
+ finalPos = 0;
+ }
+ catch(exception& e) {
+ cout << "Standard Error: " << e.what() << " has occurred in the Progress class Function Progress. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
+ exit(1);
+ }
+ catch(...) {
+ cout << "An unknown error has occurred in the Progress class function Progress. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
+ exit(1);
+ }
+}
+
/***********************************************************************/
Progress::Progress(string job, int end){
try {
cout << "********************#****#****#****#****#****#****#****#****#****#****#\n";
- cout << job << marker;
+ cout << setw(20) << left << job << setw(1) << marker;
cout.flush();
-
+
nTicks = 0;
finalPos = end;
}
/***********************************************************************/
+void Progress::newLine(string job, int end){
+ try {
+ cout << endl;
+ cout << setw(20) << left << job << setw(1) << marker;
+ cout.flush();
+
+ nTicks = 0;
+ finalPos = end;
+ }
+ catch(exception& e) {
+ cout << "Standard Error: " << e.what() << " has occurred in the Progress class Function newline. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
+ exit(1);
+ }
+ catch(...) {
+ cout << "An unknown error has occurred in the Progress class function newline. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
+ exit(1);
+ }
+}
+
+/***********************************************************************/
+
void Progress::update(const int currentPos){
try {
int ratio = int(totalTicks * (float)currentPos / finalPos);
class Progress {
public:
+ Progress();
Progress(string, int);
void update(int);
+ void newLine(string, int);
void finish();
private:
--- /dev/null
+/*
+ * slibshuff.cpp
+ * Mothur
+ *
+ * Created by Pat Schloss on 4/8/09.
+ * Copyright 2009 Patrick D. Schloss. All rights reserved.
+ *
+ */
+
+#include "slibshuff.h"
+
+/***********************************************************************/
+
+SLibshuff::SLibshuff(FullMatrix* D, int it, float co) : Libshuff(D, it, 0, co){}
+
+/***********************************************************************/
+
+float SLibshuff::evaluatePair(int i, int j){
+ return sCalculate(i,j);
+}
+
+/***********************************************************************/
+
+vector<vector<double> > SLibshuff::evaluateAll(){
+ try{
+ savedMins.resize(numGroups);
+
+ vector<vector<double> > dCXYValues(numGroups);
+
+ for(int i=0;i<numGroups;i++){
+ dCXYValues[i].resize(numGroups);
+ savedMins[i].resize(numGroups);
+ for(int j=0;j<numGroups;j++){
+ if(i!=j){
+ dCXYValues[i][j] = sCalculate(i,j);
+ savedMins[i][j] = minXY;
+ }
+
+ if(savedMins[i][i].size() == 0){
+ savedMins[i][i] = minX;
+ }
+
+ }
+ }
+ return dCXYValues;
+ }
+ catch(exception& e) {
+ cout << "Standard Error: " << e.what() << " has occurred in the SLibshuff class Function evaluateAll. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
+ exit(1);
+ }
+ catch(...) {
+ cout << "An unknown error has occurred in the SLibshuff class function evaluateAll. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
+ exit(1);
+ }
+
+}
+
+/***********************************************************************/
+
+double SLibshuff::sCalculate(int x, int y){
+ try{
+ minX = getMinX(x);
+ minXY = getMinXY(x,y);
+
+ sort(minX.begin(), minX.end());
+ sort(minXY.begin(), minXY.end());
+
+ double sum = 0.0,t=0.0;
+ int ix=0,iy=0;
+ while( (ix < groupSizes[x]) && (iy < groupSizes[x]) ) {
+ double h = (ix-iy)/double(groupSizes[x]);
+
+ if(minX[ix] < minXY[iy]) {
+ sum += (minX[ix] - t)*h*h;
+ t = minX[ix++];
+ }
+ else {
+ sum += (minXY[iy] - t)*h*h;
+ t = minXY[iy++];
+ }
+
+ }
+
+ if(ix < groupSizes[x]) {
+
+ while(ix < groupSizes[x]) {
+ double h = (ix-iy)/double(groupSizes[x]);
+ sum += (minX[ix] - t)*h*h;
+ t = minX[ix++];
+ }
+
+ }
+ else {
+
+ while(iy < groupSizes[x]) {
+ double h = (ix-iy)/double(groupSizes[x]);
+ sum += (minXY[iy] - t)*h*h;
+ t = minXY[iy++];
+ }
+
+ }
+
+ return sum;
+ }
+ catch(exception& e) {
+ cout << "Standard Error: " << e.what() << " has occurred in the SLibshuff class Function sCalculate. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
+ exit(1);
+ }
+ catch(...) {
+ cout << "An unknown error has occurred in the SLibshuff class function sCalculate. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
+ exit(1);
+ }
+
+}
+
+/***********************************************************************/
--- /dev/null
+#ifndef SLIBSHUFF
+#define SLIBSHUFF
+
+/*
+ * slibshuff.h
+ * Mothur
+ *
+ * Created by Pat Schloss on 4/8/09.
+ * Copyright 2009 Patrick D. Schloss. All rights reserved.
+ *
+ */
+
+#include "fullmatrix.h"
+#include "libshuff.h"
+
+class SLibshuff : public Libshuff {
+
+public:
+ SLibshuff(FullMatrix*, int, float);
+ vector<vector<double> > evaluateAll();
+ float evaluatePair(int, int);
+
+private:
+ double sCalculate(int, int);
+};
+
+#endif
ValidParameters::~ValidParameters() {}
-
/***********************************************************************/
bool ValidParameters::isValidParameter(string parameter, string command, string value) {
try {
try {
//{"parameter1","parameter2",...,"last parameter"};
- string readdistArray[] = {"phylip","column","name","cutoff","precision","group"};
+ string readdistArray[] = {"phylip","column", "name","cutoff","precision", "group"};
+
commandParameters["read.dist"] = addParameters(readdistArray, sizeof(readdistArray)/sizeof(string));
string readotuArray[] = {"list","order","shared", "line", "label","group","sabund", "rabund"};
projects / mothur.git / commitdiff