5 * Created by Sarah Westcott on 11/18/08.
6 * Copyright 2008 Schloss Lab UMASS Amherst. All rights reserved.
10 #include "raredisplay.h"
12 /***********************************************************************/
14 void RareDisplay::init(string label){
19 m->errorOut(e, "RareDisplay", "init");
24 /***********************************************************************/
26 void RareDisplay::update(SAbundVector* rank){
28 int newNSeqs = rank->getNumSeqs();
29 vector<double> data = estimate->getValues(rank);
33 double oldS = var[index] * ( nIters - 2 );
34 double delta = data[0] - results[index];
35 double newMean = results[index] + delta / nIters;
36 double newS = oldS + delta * ( data[0] - newMean );
37 double newVar = newS / ( nIters - 1 );
39 seqs[index] = newNSeqs;
40 results[index] = newMean;
46 seqs.push_back(newNSeqs);
47 results.push_back(data[0]);
52 m->errorOut(e, "RareDisplay", "update");
57 /***********************************************************************/
58 void RareDisplay::update(vector<SharedRAbundVector*> shared, int numSeqs, int numGroupComb) {
60 vector<double> data = estimate->getValues(shared);
61 double newNSeqs = data[0];
65 double oldS = var[index] * ( nIters - 2 );
66 double delta = data[0] - results[index];
67 double newMean = results[index] + delta / nIters;
68 double newS = oldS + delta * ( data[0] - newMean );
69 double newVar = newS / ( nIters - 1 );
70 seqs[index] = newNSeqs;
71 results[index] = newMean;
78 seqs.push_back(newNSeqs);
79 results.push_back(data[0]);
85 m->errorOut(e, "RareDisplay", "update");
90 /***********************************************************************/
92 void RareDisplay::reset(){
98 m->errorOut(e, "RareDisplay", "reset");
103 /***********************************************************************/
105 void RareDisplay::close(){
108 output->initFile(label);
110 for (int i = 0; i < seqs.size(); i++) {
112 vector<double> data(3,0);
113 double variance = var[i];
115 data[0] = results[i];
117 double ci = 1.96 * pow(variance, 0.5);
118 data[1] = data[0] - ci;
119 data[2] = data[0] + ci;
121 output->output(seqs[i], data);
133 catch(exception& e) {
134 m->errorOut(e, "RareDisplay", "close");
138 /***********************************************************************/
140 void RareDisplay::inputTempFiles(string filename){
143 m->openInputFile(filename, in);
146 in >> thisIters; m->gobble(in);
148 for (int i = 0; i < seqs.size(); i++) {
149 double tempresult, tempvar;
150 in >> tempresult >> tempvar; m->gobble(in);
152 //find weighted result
153 results[i] = ((nIters * results[i]) + (thisIters * tempresult)) / (float)(nIters + thisIters);
155 var[i] = ((nIters * var[i]) + (thisIters * tempvar)) / (float)(nIters + thisIters);
160 catch(exception& e) {
161 m->errorOut(e, "RareDisplay", "inputTempFiles");
166 /***********************************************************************/
168 void RareDisplay::outputTempFiles(string filename){
171 m->openOutputFile(filename, out);
173 out << nIters << endl;
175 for (int i = 0; i < seqs.size(); i++) {
176 out << results[i] << '\t' << var[i] << endl;
181 catch(exception& e) {
182 m->errorOut(e, "RareDisplay", "outputTempFiles");
188 /***********************************************************************/