Permalink
Cannot retrieve contributors at this time
Fetching contributors…
Cannot retrieve contributors at this time
| /* PDAP:PDTREE package for Mesquite copyright 2001-2008 P. Midford & W. MaddisonPDAP:PDTREE is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY.The web site for PDAP:PDTREE is http://mesquiteproject.org/pdap_mesquite/This source code and its compiled class files are free and modifiable under the terms of GNU Lesser General Public License. (http://www.gnu.org/copyleft/lesser.html) */package mesquite.pdap.SPOriDiagnostics;/*~~ */import JSci.maths.statistics.TDistribution;import mesquite.lib.*;import mesquite.pdap.lib.*;public class SPOriDiagnostics extends PDAP2CTStatPakSource { public boolean startJob(String arguments, Object condition, boolean hiredByName) { return true; } public PDAP2CTStatPak getStatPak(){ return new oriStatPak(); } public String getName(){ return "Regression through Origin StatPak Source"; } /*.................................................................................................................*/ public String getAuthors() { return "Peter E. Midford, Ted Garland Jr., and Wayne P. Maddison" ; } /*.................................................................................................................*/ /** returns module version */ public String getVersion() { return "1.12"; } /*.................................................................................................................*/ /** returns true when a module is prerelease */ public boolean isPrerelease() { return false; } /*.................................................................................................................*/ public boolean isSubstantive() { return false; // the statPak is substantive, but I don't think this module is... } /*.................................................................................................................*/ public String getExplanation() { return "Calculates summary statistics and regression through the origin; used for contrast vs. contrast plots"; }}class oriStatPak extends PDAP2CTStatPak{ private double r2 = MesquiteDouble.unassigned; //r squared value for least squares regression private double ts = MesquiteDouble.unassigned; //T statistic private double f = MesquiteDouble.unassigned; //F statistic private int conypos = MesquiteInteger.unassigned; //count of positive contrasts private int cony0 = 0; //count of negative contrasts private int conyneg = 0; //count of contrasts == 0 private double pvalue = MesquiteDouble.unassigned; //significance level of regression public oriStatPak() { super(); } /** positivize by flipping the sign of corresponding elements of data1 and data2 if the data1 element is negative */ private void positivize(NoPolyTree tree, double[] data1,double [] data2,int node) { for (int daughter = tree.firstDaughterOfNode(node); tree.nodeExists(daughter); daughter = tree.nextSisterOfNode(daughter)) positivize(tree, data1, data2, daughter); if (MesquiteDouble.isCombinable(data1[node]) && MesquiteDouble.isCombinable(data2[node])) { if (data1[node] < 0) { data1[node] *= -1; data2[node] *= -1; } } } private void calcBounds(Tree tree, double [] data0, double [] data1, int node) { if (tree.nodeIsInternal(node)) for (int daughter = tree.firstDaughterOfNode(node); tree.nodeExists(daughter); daughter = tree.nextSisterOfNode(daughter)) { calcBounds(tree, data0, data1, daughter); } if (MesquiteDouble.isCombinable(data0[node]) && MesquiteDouble.isCombinable(data1[node])) { if (data0[node] > max[0]) max[0] = data0[node]; if (data1[node] > max[1]) max[1] = data1[node]; if (data0[node] < min[0]) min[0] = data0[node]; if (data1[node] < min[1]) min[1] = data1[node]; } } private void calcVars(Tree tree, double[] contrasts1, double contrasts2 []) { number = 0; variance[0] = 0.0; variance[1] = 0.0; cov = 0.0; Debugg.println("(before varianceSums) number = " + number); varianceSums(tree,contrasts1,contrasts2,tree.getRoot()); Debugg.println("(after varianceSums) number = " + number); variance[0] /= getNumber()-1; variance[1] /= getNumber()-1; cov /= getNumber()-1; } /** This calculates the sums for getVars */ private void varianceSums(Tree tree, double[] contrasts1, double[] contrasts2, int node) { for (int daughter = tree.firstDaughterOfNode(node); tree.nodeExists(daughter); daughter = tree.nextSisterOfNode(daughter)) varianceSums(tree, contrasts1, contrasts2, daughter); Debugg.println("(in varianceSums; outside if) Node = " + node + "; c[1] = " + contrasts1[node] + " ; c[2] = " + contrasts2[node]); if (tree.nodeIsInternal(node) && MesquiteDouble.isCombinable(contrasts1[node]) && MesquiteDouble.isCombinable(contrasts2[node])) { number++; variance[0] += contrasts1[node]*contrasts1[node]; variance[1] += contrasts2[node]*contrasts2[node]; cov += contrasts1[node]*contrasts2[node]; } } /** @param tree the tree to traverse @param contrasts1 the contrasts in variable1 (horizontal), only used for validity checking @param contrasts2 the contrasts in variable2 (vertical), the values actually counted @param node controls recursion in the tree This method counts the number of contrasts that are greater than, less than, or equal to zero and updates the appropriate fields of the statpak. The first contrast array is just used to avoid counting points that are only defined in one variable. */ //PEM 16-Nov-2003 private void countContrasts(Tree tree, double[] contrasts1, double[] contrasts2, int node) { for (int daughter = tree.firstDaughterOfNode(node); tree.nodeExists(daughter); daughter = tree.nextSisterOfNode(daughter)) countContrasts(tree, contrasts1, contrasts2, daughter); if (tree.nodeIsInternal(node) && MesquiteDouble.isCombinable(contrasts1[node]) && MesquiteDouble.isCombinable(contrasts2[node])) { if (contrasts2[node]>0) conypos++; else if (contrasts2[node]<0) conyneg++; else cony0++; } } private double[] contrasts1, contrasts2; /** * This needs to do all the calculations; only a df reduction * count of internal nodes and ZLBs have been done at this point */ public boolean doCalculations(NumberArray xValues, NumberArray yValues, MesquiteChart chart){ min[0] = MesquiteDouble.infinite; min[1] = MesquiteDouble.infinite; max[0] = MesquiteDouble.negInfinite; max[1] = MesquiteDouble.negInfinite; resetChartBounds(chart); // this seems to be necessary for the charts to scale properly final BivariateContrastCalculator calculator = new BivariateContrastCalculator(tree, obs1, obs2); contrasts1 = calculator.contrastCalculation(); Debugg.println("After contrast calc, c1[2] = " + contrasts1[2]); contrasts2 = calculator.getContrasts2(); final int root = tree.getRoot(); positivize(tree,contrasts1,contrasts2,root); calcBounds(tree,contrasts1,contrasts2,root); Debugg.print("After calcBounds, c1[2] = " + contrasts1[2]); Debugg.println("(before calc vars) Nm1 = " + getNumberMinus1()); calcVars(tree,contrasts1,contrasts2); Debugg.println("(after calc vars) Nm1 = " + getNumberMinus1()); Debugg.print("After calcVars, c1[2] = " + contrasts1[2]); if ((variance[0] != 0) && (variance[1] != 0)) corr = cov/Math.sqrt(variance[0]*variance[1]); else corr = 0; if (variance[0] != 0) { lsregok = true; ls = cov/variance[0]; lsy = 0.0; } else lsregok = false; // reset the contrast counts conypos = 0; conyneg = 0; cony0 = 0; countContrasts(tree,contrasts1,contrasts2,root); // need xValues for filtering finishStatOrg(); tStats = new TDistribution(df); if (Math.abs(corr) <1 && df > 0) if (ts<0) pvalue= tStats.cumulative(ts); else pvalue = 1-tStats.cumulative(ts); else pvalue = MesquiteDouble.unassigned; return false; // returning true seems to cause infinite recursion in some cases } // This method finishes the calculation of the three regression models (ordinary least squares, major axis // and reduced major axis) displayed by screen 9, and now available elsewhere in Mesquite. public void finishStatOrg() { if (corr == 0) lsregok = false; if (lsregok) rma = ls/corr; if (((corr > 0) && (rma < 0)) || ((corr < 0) && (rma > 0))) rma *= -1; rmay = 0.0; final double d = Math.sqrt(((variance[0]+variance[1])*(variance[0]+variance[1])) - 4*(variance[0]*variance[1]-cov*cov)); final double lambda1 = (variance[0]+variance[1]+d)/2; // Sokal and Rolf 1981 set Y1 for vertical axis, for example on 595-596. if (lambda1 == variance[1]) { if (cov >0) ma = Double.POSITIVE_INFINITY; else ma = Double.NEGATIVE_INFINITY; } else ma = cov/(lambda1-variance[1]); // This is correct variance[1] is vertical axis may = 0.0; r2 = corr*corr; if (r2 != 1) { ts = corr*Math.sqrt((getNumberMinus1())/(1-r2)); f = ts*ts; } System.out.println("df = " + df); df = getNumberMinus1()-getDFReduction(); System.out.println("df = " + df); } public double getPValue() { return pvalue; } /*............................................................*/// Now that the statPak is set up, harvest the result string. public String flst() { if (!MesquiteDouble.isCombinable(ls)){ return ""; } StringBuffer results = new StringBuffer(300); results.append("Summary statistics for regression [of contrasts] through the origin.\n\n"); results.append("Number of contrasts: " + getNumber() + "\n"); results.append("Minimum of X,Y coordinates: " + min[0] + ", " + min[1] + "\n"); results.append("Maximum of X,Y coordinates: " + max[0] + ", " + max[1] + "\n"); results.append("Mean of X,Y coordinates: zero, zero\n"); results.append("Variance of X,Y coordinates: " + variance[0] + ", " + variance[1] + "\n"); results.append("Covariance: " + cov + "\n"); results.append("Pearson Product-Moment Correlation Coefficient: " + corr + "\n"); results.append("Reduced Major Axis."); if (lsregok) { results.append(" Slope: " + rma + "\n"); results.append(" Y Intercept: zero\n"); } else { results.append(" Slope: undefined \n"); results.append(" Y Intercept: undefined \n"); } results.append("Major Axis. \n"); results.append(" Slope: " + ma + "\n"); results.append(" Y Intercept: zero\n"); results.append("Least Squares Regression.\n"); if (lsregok) { results.append(" Slope: " + ls + "\n"); results.append(" Intercept: zero\n"); results.append(" R Squared: " + r2 + "\n"); } else { results.append(" Slope: undefined\n"); results.append(" Intercept: undefined\n"); results.append(" R Squared: undefined\n"); } if (Math.abs(corr) <1) { results.append("t: "+ ts + "\n"); results.append("F: " + f + "\n"); results.append("d.f.: " + df + "\n"); if (df > 0) { double tail1; if (ts<0) tail1= tStats.cumulative(ts); else tail1 = 1-tStats.cumulative(ts); results.append("p-value:\n"); results.append(" 2-tailed: " + (2.0*tail1) + "\n"); results.append(" 1-tailed: " + tail1 + "\n"); } } else { results.append("t: undefined\n"); results.append("F: undefined\n"); } results.append("Number of Y contrasts > 0 : " + conypos + "\n"); results.append("Number of Y contrasts < 0 : " + conyneg + "\n"); results.append("Number of Y contrasts = 0 : " + cony0 + "\n"); double signPvalue = SignTest.signTest(conypos,conyneg); if (signPvalue > 0.5) results.append("Sign test 2-tailed : 1.0\n"); else results.append("Sign test 2-tailed : " + 2.0*signPvalue + "\n"); results.append(" 1-tailed : " + signPvalue + "\n\n\n"); return results.toString(); } // This method returns text to display in the legend window. public String getLegendText() { if (!MesquiteDouble.isCombinable(ts)){ return ""; } StringBuffer results = new StringBuffer(200); results.append("Number of contrasts: " + getNumber() + "\n"); results.append("Pearson Product-Moment Correlation Coefficient: " + corr + "\n"); if (Math.abs(corr) <1) { if (df > 0) results.append("Two tailed p-value: " + (2.0*(1-tStats.cumulative(Math.abs(ts)))) + "\n"); } return results.toString(); } /** * Support for Mesquite 2.5+ 'step through trees' */ public MesquiteNumber[] getWritableResults() { if (writableResults == null){ writableResults = new MesquiteNumber[26]; writableResults[0] = new MesquiteNumber(); writableResults[0].setName("Number of contrasts"); writableResults[1] = new MesquiteNumber(); writableResults[1].setName("Min X"); writableResults[2] = new MesquiteNumber(); writableResults[2].setName("Min Y"); writableResults[3] = new MesquiteNumber(); writableResults[3].setName("Max X"); writableResults[4] = new MesquiteNumber(); writableResults[4].setName("Max Y"); writableResults[5] = new MesquiteNumber(); writableResults[5].setName("Variance X"); writableResults[6] = new MesquiteNumber(); writableResults[6].setName("Variance Y" ); writableResults[7] = new MesquiteNumber(); writableResults[7].setName("Covariance"); writableResults[8] = new MesquiteNumber(); writableResults[8].setName("Pearson Product-Moment Correlation Coefficient"); writableResults[9] = new MesquiteNumber(); writableResults[9].setName("RMA Slope"); writableResults[10] = new MesquiteNumber(); writableResults[10].setName("RMA Intercept"); writableResults[11] = new MesquiteNumber(); writableResults[11].setName("MA Slope"); writableResults[12] = new MesquiteNumber(); writableResults[12].setName("MA Intercept"); writableResults[13] = new MesquiteNumber(); writableResults[13].setName("OLS Slope"); writableResults[14] = new MesquiteNumber(); writableResults[14].setName("OLS Intercept"); writableResults[15] = new MesquiteNumber(); writableResults[15].setName("r-squared"); writableResults[16] = new MesquiteNumber(); writableResults[16].setName("t"); writableResults[17] = new MesquiteNumber(); writableResults[17].setName("F"); writableResults[18] = new MesquiteNumber(); writableResults[18].setName("df"); writableResults[19] = new MesquiteNumber(); writableResults[19].setName("p (2-tailed)"); writableResults[20] = new MesquiteNumber(); writableResults[20].setName("p (1-tailed)"); writableResults[21] = new MesquiteNumber(); writableResults[21].setName("y contrasts > 0"); writableResults[22] = new MesquiteNumber(); writableResults[22].setName("y contrasts < 0"); writableResults[23] = new MesquiteNumber(); writableResults[23].setName("y contrasts = 0"); writableResults[24] = new MesquiteNumber(); writableResults[24].setName("sign test (2-tailed)"); writableResults[25] = new MesquiteNumber(); writableResults[25].setName("sign test (1-tailed)"); } if (MesquiteDouble.isCombinable(corr)){ writableResults[0].setValue(getNumber()); writableResults[1].setValue(min[0]); writableResults[2].setValue(min[1]); writableResults[3].setValue(max[0]); writableResults[4].setValue(max[1]); writableResults[5].setValue(variance[0]); writableResults[6].setValue(variance[1]); writableResults[7].setValue(cov); writableResults[8].setValue(corr); writableResults[9].setValue(rma); writableResults[10].setValue(rmay); writableResults[11].setValue(ma); writableResults[12].setValue(may); writableResults[13].setValue(ls); writableResults[14].setValue(lsy); writableResults[15].setValue(r2); writableResults[16].setValue(ts); writableResults[17].setValue(f); writableResults[18].setValue(df); if (df > 0) { double tail1; if (ts<0) tail1= tStats.cumulative(ts); else tail1 = 1-tStats.cumulative(ts); writableResults[19].setValue(2*tail1); writableResults[20].setValue(tail1); } else{ writableResults[19].setToUnassigned(); writableResults[20].setToUnassigned(); } writableResults[21].setValue(conypos); writableResults[22].setValue(conyneg); writableResults[23].setValue(cony0); final double signPvalue = SignTest.signTest(conypos,conyneg); if (signPvalue > 0.5){ writableResults[24].setValue(1.0); writableResults[25].setValue(signPvalue); } else{ writableResults[24].setValue(2.0*signPvalue); writableResults[25].setValue(signPvalue); } } else for (int i=0;i<writableResults.length;i++) writableResults[i].setToUnassigned(); return writableResults; } public String CIString(String lineEnding, String tableDelimiter, boolean convertSpaces) { return ""; }} |