1 #include "trialswap2.h"
4 //The sum_of_squares, havel_hakimi and calc_c_score algorithms have been adapted from I. Miklos and J. Podani. 2004. Randomization of presence-absence matrices: comments and new algorithms. Ecology 85:86-92.
7 /**************************************************************************************************
8 int TrialSwap2::intrand(int n){
12 z = (double)random() * (double)n / (double)RAND_MAX;
17 return((int)floor(z));
20 m->errorOut(e, "TrialSwap2", "intrand");
24 /**************************************************************************************************/
25 /* completely random matrix, all column and row totals are variable, matrix size is the same
29 /**************************************************************************************************/
30 int TrialSwap2::sim1(vector<vector<int> > &co_matrix){
33 vector<vector<int> > tmpmatrix;
34 int nrows = co_matrix.size();
35 int ncols = co_matrix[0].size();
38 // for(i=0;i<nrows;i++)
43 //cout << "building matrix" << endl;
44 for(int i=0;i<nrows;i++){
45 if (m->control_pressed) { break; }
47 for(int j=0;j<ncols;j++){
48 double randNum = rand() / double(RAND_MAX);
49 //cout << randNum << endl;
57 tmpmatrix.push_back(randRow);
61 co_matrix = tmpmatrix;
66 m->errorOut(e, "TrialSwap2", "sim1");
70 /**************************************************************************************************/
72 *row sums fixed, columns equiprobable
74 void TrialSwap2::sim2(vector<vector<int> > &co_matrix)
78 for(int i=0;i<co_matrix.size();i++)
80 if (m->control_pressed) { break; }
81 random_shuffle( co_matrix[i].begin(), co_matrix[i].end() );
85 m->errorOut(e, "TrialSwap2", "sim2");
89 /**************************************************************************************************/
90 int TrialSwap2::sim2plus(vector<int> rowtotal, vector<vector<int> > &co_matrix)
93 int nrows = co_matrix.size();
94 int ncols = co_matrix[0].size();
95 double cellprob = 1.0/ncols;
96 vector<double> cellprobvec;
98 vector<vector<int> > tmpmatrix;
103 for(int i=0; i<ncols; i++)
105 if (m->control_pressed) { return 0; }
106 cellprobvec.push_back(start + cellprob);
107 start = cellprobvec[i];
110 for(int i=0; i<nrows; i++)
112 tmprow.assign(ncols, 0);
114 while( accumulate( tmprow.begin(), tmprow.end(), 0 ) < rowtotal[i])
116 if (m->control_pressed) { return 0; }
117 double randNum = rand() / double(RAND_MAX);
118 //cout << randNum << endl;
119 if(randNum <= cellprobvec[0])
124 for(int j=1;j<ncols;j++)
126 //cout << range[j] << endl;
127 if(randNum <= cellprobvec[j] && randNum > cellprobvec[j-1] && tmprow[j] != 1)
133 tmpmatrix.push_back(tmprow);
136 co_matrix = tmpmatrix;
142 catch(exception& e) {
143 m->errorOut(e, "TrialSwap2", "sim2plus");
147 /**************************************************************************************************/
149 * same as sim2 but using initmatrix which is the initial co-occurrence matrix before transposition
150 * may have to be changed depending on what matrix 'seed' is used. One way to use is to transpose
151 * every null matrix before using an index and use the random matrix as a seed for the next null.
153 /**************************************************************************************************/
154 void TrialSwap2::sim3(vector<vector<int> > &initmatrix)
157 for(int i=0;i<initmatrix.size();i++)
159 if (m->control_pressed) { break; }
160 random_shuffle( initmatrix[i].begin(), initmatrix[i].end() );
164 catch(exception& e) {
165 m->errorOut(e, "TrialSwap2", "sim3");
169 /**************************************************************************************************/
175 /**************************************************************************************************/
176 int TrialSwap2::sim4(vector<int> columntotal, vector<int> rowtotal, vector<vector<int> > &co_matrix)
179 vector<double> colProb;
180 vector<int> tmprow;//(ncols, 7);
181 vector<vector<int> > tmpmatrix;
182 vector<double> range;
183 vector<double> randNums;
184 int ncols = columntotal.size();
185 int nrows = rowtotal.size();
188 double colSum = accumulate( columntotal.begin(), columntotal.end(), 0 );
189 //cout << "col sum: " << colSum << endl;
190 for(int i=0;i<ncols;i++)
192 if (m->control_pressed) { return 0; }
193 colProb.push_back(columntotal[i]/colSum);
198 for(int i=0;i<ncols;i++)
200 if (m->control_pressed) { return 0; }
201 range.push_back(start + colProb[i]);
205 for(int i=0;i<nrows;i++)
207 tmprow.assign(ncols, 0);
208 if (m->control_pressed) { return 0; }
210 while ( accumulate( tmprow.begin(), tmprow.end(), 0 ) < rowtotal[i])
212 if (m->control_pressed) { return 0; }
214 double randNum = rand() / double(RAND_MAX);
215 if(randNum <= range[0])
220 for(int j=1;j<ncols;j++)
222 if(randNum <= range[j] && randNum > range[j-1] && tmprow[j] != 1)
229 tmpmatrix.push_back(tmprow);
233 co_matrix = tmpmatrix;
237 catch(exception& e) {
238 m->errorOut(e, "TrialSwap2", "sim4");
242 /**************************************************************************************************/
244 * inverse of sim4, MUST BE TRANSPOSED BEFORE CO-OCCURRENCE ANALYSIS
248 /**************************************************************************************************/
249 int TrialSwap2::sim5(vector<int> initcolumntotal,vector<int> initrowtotal, vector<vector<int> > &initmatrix)
252 vector<double> colProb;
253 vector<int> tmprow;//(ncols, 7);
254 vector<vector<int> > tmpmatrix;
255 vector<double> range;
256 vector<double> randNums;
257 int ncols = initcolumntotal.size();
258 int nrows = initrowtotal.size();
262 double colSum = accumulate( initcolumntotal.begin(), initcolumntotal.end(), 0 );
263 //cout << "col sum: " << colSum << endl;
264 for(int i=0;i<ncols;i++)
266 if (m->control_pressed) { return 0; }
267 colProb.push_back(initcolumntotal[i]/colSum);
272 for(int i=0;i<ncols;i++)
274 if (m->control_pressed) { return 0; }
275 range.push_back(start + colProb[i]);
279 for(int i=0;i<nrows;i++)
281 tmprow.assign(ncols, 0);
282 if (m->control_pressed) { return 0; }
284 while ( accumulate( tmprow.begin(), tmprow.end(), 0 ) < initrowtotal[i])
286 if (m->control_pressed) { return 0; }
288 double randNum = rand() / double(RAND_MAX);
289 if(randNum <= range[0])
294 for(int j=1;j<ncols;j++)
296 if(randNum <= range[j] && randNum > range[j-1] && tmprow[j] != 1)
303 tmpmatrix.push_back(tmprow);
307 initmatrix = tmpmatrix;
310 catch(exception& e) {
311 m->errorOut(e, "TrialSwap2", "sim5");
315 /**************************************************************************************************/
321 /**************************************************************************************************/
322 int TrialSwap2::sim6(vector<int> columntotal, vector<vector<int> > &co_matrix)
325 vector<vector<int> > tmpmatrix;
326 vector<double> colProb;
328 vector<double> range;
329 int ncols = columntotal.size();
330 int nrows = co_matrix.size();
332 int colSum = accumulate( columntotal.begin(), columntotal.end(), 0 );
334 for(int i=0;i<ncols;i++)
336 if (m->control_pressed) { return 0; }
337 colProb.push_back(columntotal[i]/double (colSum));
342 for(int i=0;i<ncols;i++)
344 if (m->control_pressed) { return 0; }
345 range.push_back(start + colProb[i]);
349 for(int i=0;i<nrows;i++)
351 if (m->control_pressed) { return 0; }
352 tmprow.assign(ncols, 0);
354 tmprowtotal = (rand() / double (RAND_MAX)) * 10;
355 while ( tmprowtotal > ncols) {
356 if (m->control_pressed) { return 0; }
357 tmprowtotal = (rand() / double (RAND_MAX)) * 10;
359 //cout << tmprowtotal << endl;
360 //cout << accumulate( tmprow.begin(), tmprow.end(), 0 ) << endl;
362 while ( accumulate( tmprow.begin(), tmprow.end(), 0 ) < tmprowtotal)
364 if (m->control_pressed) { return 0; }
365 double randNum = rand() / double(RAND_MAX);
366 //cout << randNum << endl;
367 if(randNum <= range[0])
372 for(int j=1;j<ncols;j++)
374 //cout << range[j] << endl;
375 if(randNum <= range[j] && randNum > range[j-1] && tmprow[j] != 1)
385 tmpmatrix.push_back(tmprow);
389 co_matrix = tmpmatrix;
394 catch(exception& e) {
395 m->errorOut(e, "TrialSwap2", "sim6");
399 /**************************************************************************************************/
401 * MUST BE TRANSPOSED BEFORE CO-OCCURRENCE ANALYSIS
405 /**************************************************************************************************/
406 int TrialSwap2::sim7(vector<int> initrowtotal, vector<vector<int> > &co_matrix)
409 vector<vector<double> > probmatrix;
410 vector<vector<int> > tmpmatrix;
411 vector<double> colProb;
412 vector<double> probrow;
414 vector<double> range;
416 int ncols = co_matrix[0].size(); int nrows = co_matrix.size();
418 tmpmatrix.assign(nrows, vector<int>(ncols, 0.));
420 int rowsum = accumulate( initrowtotal.begin(), initrowtotal.end(), 0 );
423 //cout << nc << endl;
425 //assign null matrix based on probabilities
427 double start = 0.0; // don't reset start -- probs should be from 0-1 thoughout the entire matrix
429 for(int i=0;i<nrows;i++)
431 if (m->control_pressed) { return 0; }
432 //cout << initrowtotal[i]/double(nc) << endl;
433 double cellprob = initrowtotal[i]/double(nc);
434 //cout << cellprob << endl;
435 for(int j=0;j<ncols;j++)
438 probrow.push_back(start + cellprob);
439 //cout << probrow[j] << endl;
440 //cout << start << endl;
441 start = start + cellprob;
443 probmatrix.push_back(probrow);
448 //while(tmprowsum < rowsum)
449 //for(int k=0;k<rowsum;k++)
453 if (m->control_pressed) { return 0; }
456 //tmprowsum = accumulate( tmprowtotal.begin(), tmprowtotal.end(), 0 );
457 double randNum = rand() / double(RAND_MAX);
458 //cout << randNum << "+" << endl;
459 //special case for the first entry
460 if(randNum <= probmatrix[0][0] && tmpmatrix[0][0] != 1)
469 for(int i=0;i<nrows;i++)
471 if (m->control_pressed) { return 0; }
472 for(int j=0;j<ncols;j++)
474 //cout << probmatrix[i][j] << endl;
475 if(randNum <= probmatrix[i][j] && randNum > probmatrix[i][j-1] && tmpmatrix[i][j] != 1)
490 co_matrix = tmpmatrix;
492 //build probibility matrix
493 /* for(int i=0;i<nrows;i++)
495 for(int j=0;j<ncols;j++)
497 probrow.push_back(rowtotal[i]/nc);
499 probmatrix.pushback(probrow);
504 /* int colSum = accumulate( initcolumntotal.begin(), initcolumntotal.end(), 0 );
506 for(int i=0;i<ncols;i++)
508 colProb.push_back(initcolumntotal[i]/double (colSum));
513 for(int i=0;i<ncols;i++)
515 range.push_back(start + colProb[i]);
519 for(int i=0;i<nrows;i++)
521 tmprow.assign(ncols, 0);
523 tmprowtotal = (rand() / double (RAND_MAX)) * 10;
524 while ( tmprowtotal > ncols)
525 tmprowtotal = (rand() / double (RAND_MAX)) * 10;
526 //cout << tmprowtotal << endl;
527 //cout << accumulate( tmprow.begin(), tmprow.end(), 0 ) << endl;
529 while ( accumulate( tmprow.begin(), tmprow.end(), 0 ) < tmprowtotal)
531 double randNum = rand() / double(RAND_MAX);
532 //cout << randNum << endl;
533 if(randNum <= range[0])
538 for(int j=1;j<ncols;j++)
540 //cout << range[j] << endl;
541 if(randNum <= range[j] && randNum > range[j-1] && tmprow[j] != 1)
548 tmpmatrix.push_back(tmprow);
552 initmatrix = tmpmatrix;
555 catch(exception& e) {
556 m->errorOut(e, "TrialSwap2", "sim7");
560 /**************************************************************************************************/
566 /**************************************************************************************************/
567 int TrialSwap2::sim8(vector<int> columntotal, vector<int> rowtotal, vector<vector<int> > &co_matrix)
572 int ncols = columntotal.size(); int nrows = rowtotal.size();
573 double probarray[nrows * ncols];
578 //double colSum = accumulate( columntotal.begin(), columntotal.end(), 0 );
579 double rowSum = accumulate( rowtotal.begin(), rowtotal.end(), 0 );
581 if (m->control_pressed) { return 0; }
583 //cout << "rowsum: " << rowSum << endl;
587 //create probability matrix with each site being between 0 and 1
588 for (int i=0;i<nrows;i++) {
589 if (m->control_pressed) { return 0; }
590 for (int j=0;j<ncols;j++) {
591 prob = (rowtotal[i] * columntotal[j])/(rowSum*rowSum);
593 probarray[ncols * i + j] = -1;
595 probarray[ncols * i + j] = start + prob;
596 //probmatrixrow.pushback(start + prob);
600 //cout << "prbarray" << endl;
601 //for(int i=0;i<(nrows*ncols);i++)
602 //cout << probarray[i] << " ";
605 //generate random muber between 0 and 1 and interate through probarray until found
606 while (total < grandtotal) {
607 if (m->control_pressed) { return 0; }
608 randnum = rand() / double(RAND_MAX);
609 //cout << "rand num: " << randnum << endl;
610 if((randnum <= probarray[0]) && (probarray[0] != 2) ) {
615 for(int i=1;i<(nrows*ncols);i++) {
616 if (m->control_pressed) { return 0; }
617 if((randnum <= probarray[i]) && (randnum > probarray[i-1]) && (probarray[i] != 2) ) {
626 //cout << "prbarray" << endl;
627 //for(int i=0;i<(nrows*ncols);i++)
628 //cout << probarray[i] << " ";
630 for(int i=0;i<nrows;i++) {
631 if (m->control_pressed) { return 0; }
632 for(int j=0;j<ncols;j++) {
633 if(probarray[ncols * i + j] == 2)
641 catch(exception& e) {
642 m->errorOut(e, "TrialSwap2", "sim8");
646 /**************************************************************************************************/
647 double TrialSwap2::calc_c_score (vector<vector<int> > &co_matrix,vector<int> rowtotal)
653 double normcscore = 0.0;
655 int ncols = co_matrix[0].size(); int nrows = rowtotal.size();
656 vector<vector<double> > s; s.resize(nrows);
657 for (int i = 0; i < nrows; i++) { s[i].resize(nrows,0.0); }//only fill half the matrix
660 for(int i=0;i<nrows-1;i++)
663 for(int j=i+1;j<nrows;j++)
665 if (m->control_pressed) { return 0; }
666 for(int k=0;k<ncols;k++)
668 if((co_matrix[i][k]==1)&&(co_matrix[j][k]==1)) //if both are 1s ie co-occurrence
669 s[i][j]++; //s counts co-occurrences
672 //rowtotal[i] = A, rowtotal[j] = B, ncols = P, s[i][j] = J
673 cscore += (rowtotal[i]-s[i][j])*(rowtotal[j]-s[i][j]);///(nrows*(nrows-1)/2);
674 D = (rowtotal[i]-s[i][j])*(rowtotal[j]-s[i][j]);
676 if(ncols < (rowtotal[i] + rowtotal[j]))
678 maxD = (ncols-rowtotal[i])*(ncols-rowtotal[j]);
682 maxD = rowtotal[i] * rowtotal[j];
687 normcscore += D/maxD;
693 cscore = cscore/(double)(nrows*(nrows-1)/2);
694 //cout << "normalized c score: " << normcscore/nonzeros << endl;
698 catch(exception& e) {
699 m->errorOut(e, "TrialSwap2", "calc_c_score");
703 /**************************************************************************************************/
704 int TrialSwap2::calc_checker (vector<vector<int> > &co_matrix, vector<int> rowtotal)
708 //int s[nrows][ncols];
709 int ncols = co_matrix[0].size(); int nrows = rowtotal.size();
710 vector<vector<int> > s; s.resize(nrows);
711 for (int i = 0; i < nrows; i++) { s[i].resize(nrows,0); }//only fill half the matrix
713 for(int i=0;i<nrows-1;i++)
715 for(int j=i+1;j<nrows;j++)
717 if (m->control_pressed) { return 0; }
719 for(int k=0;k<ncols;k++)
721 //cout << s[i][j] << endl;
722 //iterates through the row and counts co-occurrences. The total number of co-occurrences for each row pair is kept in matrix s at location s[i][j].
723 if((co_matrix[i][k]==1)&&(co_matrix[j][k]==1)) //if both are 1s ie co-occurrence
724 s[i][j]++; //s counts co-occurrences
727 //cout << "rowtotal: " << rowtotal[i] << endl;
728 //cout << "co-occurrences: " << s[i][j] << endl;
729 //cunits+=(rowtotal[i]-s[i][j])*(rowtotal[j]-s[i][j]);
740 catch(exception& e) {
741 m->errorOut(e, "TrialSwap2", "calc_checker");
745 /**************************************************************************************************/
746 double TrialSwap2::calc_vratio (vector<int> rowtotal, vector<int> columntotal)
749 int nrows = rowtotal.size();
750 int ncols = columntotal.size();
751 int sumCol = accumulate(columntotal.begin(), columntotal.end(), 0 );
752 // int sumRow = accumulate(rowtotal.begin(), rowtotal.end(), 0 );
754 double colAvg = (double) sumCol / (double) ncols;
755 // double rowAvg = (double) sumRow / (double) nrows;
759 // double totalRowVar = 0.0;
763 for(int i=0;i<nrows;i++)
765 if (m->control_pressed) { return 0; }
766 p = (double) rowtotal[i]/(double) ncols;
767 rowVar += p * (1.0-p);
770 for(int i=0;i<ncols;i++)
772 if (m->control_pressed) { return 0; }
773 colVar += pow(((double) columntotal[i]-colAvg),2);
776 colVar = (1.0/(double)ncols) * colVar;
778 return colVar/rowVar;
780 catch(exception& e) {
781 m->errorOut(e, "TrialSwap2", "calc_vratio");
786 /**************************************************************************************************/
787 int TrialSwap2::calc_combo (vector<vector<int> > &initmatrix)
790 int initrows = initmatrix.size();
794 for(int i=0;i<initrows;i++)
797 for(int j=i+1;j<=initrows;j++)
799 if (m->control_pressed) { return 0; }
800 if( (initmatrix[i] == initmatrix[j]))
808 //on the last iteration of a previously matched row it will add itself because it doesn't match any following rows, so that combination is counted
814 catch(exception& e) {
815 m->errorOut(e, "TrialSwap2", "calc_combo");
819 /**************************************************************************************************/
820 int TrialSwap2::swap_checkerboards (vector<vector<int> > &co_matrix)
823 int ncols = co_matrix[0].size(); int nrows = co_matrix.size();
825 i = m->getRandomIndex(nrows-1);
826 while((j = m->getRandomIndex(nrows-1) ) == i ) {;if (m->control_pressed) { return 0; }}
827 k = m->getRandomIndex(ncols-1);
828 while((l = m->getRandomIndex(ncols-1)) == k ) {;if (m->control_pressed) { return 0; }}
830 //cout << co_matrix[i][k] << " " << co_matrix[j][l] << endl;
831 //cout << co_matrix[i][l] << " " << co_matrix[j][k] << endl;
832 //cout << co_matrix[i][l] << " " << co_matrix[j][k] << endl;
833 //cout << co_matrix[i][l] << " " << co_matrix[j][k] << endl;
834 if((co_matrix[i][k]*co_matrix[j][l]==1 && co_matrix[i][l]+co_matrix[j][k]==0)||(co_matrix[i][k]+co_matrix[j][l]==0 && co_matrix[i][l]*co_matrix[j][k]==1)) //checking for checkerboard value and swap
836 co_matrix[i][k]=1-co_matrix[i][k];
837 co_matrix[i][l]=1-co_matrix[i][l];
838 co_matrix[j][k]=1-co_matrix[j][k];
839 co_matrix[j][l]=1-co_matrix[j][l];
840 //cout << "swapped!" << endl;
842 //cout << "i: " << i << " j: " << j << " k: " << " l: " << l << endl;
845 catch(exception& e) {
846 m->errorOut(e, "TrialSwap2", "swap_checkerboards");
850 /**************************************************************************************************/
851 double TrialSwap2::calc_pvalue_greaterthan (vector<double> scorevec, double initialscore)
854 int runs = scorevec.size();
856 for( int i=0;i<runs;i++)
858 if (m->control_pressed) { return 0; }
859 if(scorevec[i]>=initialscore)
862 return p/(double)runs;
864 catch(exception& e) {
865 m->errorOut(e, "TrialSwap2", "calc_pvalue_greaterthan");
869 /**************************************************************************************************/
870 double TrialSwap2::calc_pvalue_lessthan (vector<double> scorevec, double initialscore)
873 int runs = scorevec.size();
875 for( int i=0;i<runs;i++)
877 if (m->control_pressed) { return 0; }
878 if(scorevec[i]<=initialscore)
881 return p/(double)runs;
883 catch(exception& e) {
884 m->errorOut(e, "TrialSwap2", "calc_pvalue_lessthan");
888 /**************************************************************************************************/
889 double TrialSwap2::t_test (double initialscore, int runs, double nullMean, vector<double> scorevec)
896 for(int i=0;i<runs;i++)
898 if (m->control_pressed) { return 0; }
899 sum += pow((scorevec[i] - nullMean),2);
900 //cout << "scorevec[" << i << "]" << scorevec[i] << endl;
903 m->mothurOut("nullMean: " + toString(nullMean)); m->mothurOutEndLine();
905 m->mothurOut("sum: " + toString(sum)); m->mothurOutEndLine();
907 sampleSD = sqrt( (1/runs) * sum );
909 m->mothurOut("samplSD: " + toString(sampleSD)); m->mothurOutEndLine();
911 t = (nullMean - initialscore) / (sampleSD / sqrt(runs));
915 catch(exception& e) {
916 m->errorOut(e, "TrialSwap2", "t_test");
920 /**************************************************************************************************/
921 int TrialSwap2::print_matrix(vector<vector<int> > &matrix, int nrows, int ncols)
924 m->mothurOut("matrix:"); m->mothurOutEndLine();
926 for (int i = 0; i < nrows; i++)
928 if (m->control_pressed) { return 0; }
929 for (int j = 0; j < ncols; j++)
931 m->mothurOut(toString(matrix[i][j]));
933 m->mothurOutEndLine();
937 catch(exception& e) {
938 m->errorOut(e, "TrialSwap2", "print_matrix");
942 /**************************************************************************************************/
943 int TrialSwap2::transpose_matrix (vector<vector<int> > &initmatrix, vector<vector<int> > &co_matrix)//, int nrows, int nocols)
946 int ncols = initmatrix.size(); int nrows = initmatrix[0].size();
947 int tmpnrows = nrows;
948 //vector<vector<int> > tmpvec;
950 if(!co_matrix.empty())
952 for (int i=0;i<nrows;i++)
954 if (m->control_pressed) { return 0; }
955 for (int j=0;j<ncols;j++)
957 tmprow.push_back(initmatrix[j][i]);
959 /*if (accumulate( tmprow.begin(), tmprow.end(), 0 ) == 0)
964 co_matrix.push_back(tmprow);
970 catch(exception& e) {
971 m->errorOut(e, "TrialSwap2", "transpose_matrix");
975 /**************************************************************************************************/
976 int TrialSwap2::update_row_col_totals(vector<vector<int> > &co_matrix, vector<int> &rowtotal, vector<int> &columntotal)
980 //columntotal.clear();
981 //generate (rowtotal.begin(), rowtotal.end(), 0);
982 //generate (columntotal.begin(), columntotal.end(), 0);
983 int nrows = co_matrix.size();
984 int ncols = co_matrix[0].size();
985 vector<int> tmpcolumntotal; tmpcolumntotal.resize(ncols, 0);
986 vector<int> tmprowtotal; tmprowtotal.resize(nrows, 0);
990 for (int i = 0; i < nrows; i++)
992 if (m->control_pressed) { return 0; }
993 for (int j = 0; j < ncols; j++)
995 if (co_matrix[i][j] == 1)
1001 tmprowtotal[i] = rowcount;
1004 columntotal = tmpcolumntotal;
1005 rowtotal = tmprowtotal;
1006 /*cout << "rowtotal: ";
1007 for(int i = 0; i<nrows; i++) { cout << rowtotal[i]; }
1009 cout << " coltotal: ";
1010 for(int i = 0; i<ncols; i++) { cout << columntotal[i]; }
1014 catch(exception& e) {
1015 m->errorOut(e, "TrialSwap2", "update_row_col_totals");
1019 /**************************************************************************************************/