001/* 002 * BioJava development code 003 * 004 * This code may be freely distributed and modified under the 005 * terms of the GNU Lesser General Public Licence. This should 006 * be distributed with the code. If you do not have a copy, 007 * see: 008 * 009 * http://www.gnu.org/copyleft/lesser.html 010 * 011 * Copyright for this code is held jointly by the individual 012 * authors. These should be listed in @author doc comments. 013 * 014 * For more information on the BioJava project and its aims, 015 * or to join the biojava-l mailing list, visit the home page 016 * at: 017 * 018 * http://www.biojava.org/ 019 * 020 */ 021package org.biojava.nbio.survival.cox; 022 023import org.biojava.nbio.survival.cox.stats.ChiSq; 024import org.biojava.nbio.survival.kaplanmeier.figure.ExpressionFigure; 025import org.biojava.nbio.survival.kaplanmeier.figure.KaplanMeierFigure; 026 027import java.text.DecimalFormat; 028import java.util.ArrayList; 029import java.util.LinkedHashMap; 030 031//import org.biojava.nbio.survival.cox.comparators.SurvivalInfoComparator; 032 033/** 034 * Holds the results of a cox analysis where calling dump(), toString() will give an output similar to R 035 * @author Scooter Willis <willishf at gmail dot com> 036 */ 037public class CoxInfo { 038 039 private WaldTestInfo waldTestInfo = null; 040 String message = ""; 041 Integer maxIterations = null; 042 Double eps = null; 043 Double toler = null; 044 CoxMethod method; 045 private double[][] imat = null; //:the variance matrix at beta=final 046 private double[][] naive_imat = null; //the original variance matrix used in residuals 047 double[] u = new double[0]; //:score vector 048 int iterations = 0; //:actual number of iterations used 049 int flag = 0; // success flag 1000 did not converge 1 to nvar: rank of the solution 050 double logTest = 0; 051 double logTestpval = 0; 052 double loglikInit = 0; //loglik at beta=initial values, at beta=final 053 double loglikFinal = 0; 054 Double scoreLogrankTest; 055 private Double rscore = null; //robust score 056 private Double rscoreLogrankTestpvalue = null; 057 private double degreeFreedom; 058 private Double scoreLogrankTestpvalue; 059 int numSamples = 0; 060 int numEvents = 0; 061 private LinkedHashMap<String, String> metaDataFilter = null; 062 private LinkedHashMap<String, CoxCoefficient> coefficientsList = new LinkedHashMap<String, CoxCoefficient>(); 063 LinkedHashMap<Double, Double> baselineSurvivorFunction = new LinkedHashMap<Double, Double>(); 064 ArrayList<SurvivalInfo> survivalInfoList = new ArrayList<SurvivalInfo>(); 065 /** 066 * 067 */ 068 public KaplanMeierFigure kmf = null; 069 /** 070 * 071 */ 072 public ExpressionFigure ef = null; 073 074 /** 075 * 076 * @return 077 */ 078 public ArrayList<SurvivalInfo> getSurvivalInfoList() { 079 return survivalInfoList; 080 } 081 082 /** 083 * 084 * @param var 085 * @throws Exception 086 */ 087 public void setVariance(double[][] var) throws Exception { 088 089 // if (Math.abs(var[0][1] - var[1][0]) > .0000000000001) { //in the CoxCC correction looks like a precision error keeps these from being equal 10-19 090 // throw new Exception("Expecting diagonal to be equal"); 091 // } 092 imat = new double[var.length][var[0].length]; 093 for (int i = 0; i < var.length; i++) { 094 for (int j = 0; j < var[0].length; j++) { 095 imat[i][j] = var[i][j]; 096 } 097 } 098 calcSummaryValues(); 099 } 100 101 /** 102 * 103 * @return 104 */ 105 public double[][] getVariance() { 106 double[][] var = new double[imat.length][imat[0].length]; 107 for (int i = 0; i < var.length; i++) { 108 for (int j = 0; j < var[0].length; j++) { 109 var[i][j] = imat[i][j]; 110 } 111 } 112 113 return var; 114 } 115 116 /** 117 * 118 * @param var 119 * @throws Exception 120 */ 121 public void setNaiveVariance(double[][] var) throws Exception { 122// if (var[0][1] != var[1][0]) { 123// throw new Exception("Expecting diagonal to be equal"); 124// } 125 126 127 naive_imat = new double[var.length][var[0].length]; 128 for (int i = 0; i < var.length; i++) { 129 for (int j = 0; j < var[0].length; j++) { 130 naive_imat[i][j] = var[i][j]; 131 } 132 } 133 134 calcSummaryValues(); 135 } 136 137 /** 138 * 139 * @return 140 */ 141 public double[][] getNaiveVariance() { 142 double[][] var = new double[imat.length][imat[0].length]; 143 for (int i = 0; i < var.length; i++) { 144 for (int j = 0; j < var[0].length; j++) { 145 var[i][j] = naive_imat[i][j]; 146 } 147 } 148 149 return var; 150 151 } 152 153 /** 154 * 155 * @param data 156 */ 157 public void setSurvivalInfoList(ArrayList<SurvivalInfo> data) { 158 survivalInfoList = data; 159 numSamples = data.size(); 160 161 for (SurvivalInfo si : data) { 162 if (si.getStatus() == 1) { 163 numEvents++; 164 } 165 } 166 } 167 168 /** 169 * 170 * @return 171 */ 172 public double[] getWeighted() { 173 double[] weighted = new double[survivalInfoList.size()]; 174 int p = 0; 175 for (SurvivalInfo si : this.survivalInfoList) { 176 weighted[p] = si.getWeight(); 177 p++; 178 } 179 return weighted; 180 } 181 182 /** 183 * 184 * @return 185 */ 186 public double[][] getVariableResiduals() { 187 ArrayList<String> variables = new ArrayList<String>(coefficientsList.keySet()); 188 double[][] rr = new double[survivalInfoList.size()][variables.size()]; 189 int p = 0; 190 for (SurvivalInfo si : this.survivalInfoList) { 191 int i = 0; 192 for (String v : variables) { 193 rr[p][i] = si.getResidualVariable(v); 194 i++; 195 } 196 p++; 197 } 198 199 return rr; 200 } 201 202 /** 203 * 204 * @param rr 205 */ 206 public void setVariableResiduals(double[][] rr) { 207 ArrayList<String> variables = new ArrayList<String>(coefficientsList.keySet()); 208 209 int p = 0; 210 for (SurvivalInfo si : this.survivalInfoList) { 211 int i = 0; 212 for (String v : variables) { 213 si.setResidualVariable(v, rr[p][i]); 214 i++; 215 } 216 p++; 217 } 218 219 } 220 221 /** 222 * 223 * @return 224 */ 225 public int getNumberCoefficients() { 226 return coefficientsList.size(); 227 } 228 229 /** 230 * 231 * @param name 232 * @return 233 */ 234 public CoxCoefficient getCoefficient(String name) { 235 return coefficientsList.get(name); 236 } 237 238 /** 239 * 240 * @param name 241 * @param coefficient 242 */ 243 public void setCoefficient(String name, CoxCoefficient coefficient) { 244 coefficientsList.put(name, coefficient); 245 } 246 247 /** 248 * 249 * @param header 250 * @param beginLine 251 * @param beginCell 252 * @param endCell 253 * @param endLine 254 * @return 255 */ 256 public String getCoefficientText(boolean header, String beginLine, String beginCell, String endCell, String endLine) { 257 String o = ""; 258 if (header) { 259 String robust = ""; 260 if (naive_imat != null) { 261 robust = beginCell + "robust se" + endCell; 262 } 263 o = o + beginLine + beginCell + fmtpl("", 9) + endCell + beginCell + fmtpl("coef", 9) + endCell + beginCell + fmtpl("se(coef)", 9) + endCell + robust + beginCell + fmtpl("z", 9) + endCell + beginCell + fmtpl("p-value", 9) + endCell + beginCell + fmtpl("HR", 9) + endCell + beginCell + fmtpl("lower .95", 9) + endCell + beginCell + fmtpl("upper .95", 9) + endCell + endLine; 264 }//Coefficients,Coe,StdErr,HR,p-value,HR Lo 95%,HR Hi 95% 265 266 for (CoxCoefficient coe : coefficientsList.values()) { 267 String robust = ""; 268 String stderror = ""; 269 if (naive_imat != null) { 270 stderror = beginCell + fmt(coe.getRobustStdError(), 5, 9) + endCell; 271 robust = beginCell + fmt(coe.getStdError(), 5, 9) + endCell; 272 } else { 273 stderror = beginCell + fmt(coe.getStdError(), 5, 9) + endCell; 274 } 275 o = o + beginLine + beginCell + fmtpr(coe.getName(), 9) + endCell + beginCell + fmt(coe.getCoeff(), 5, 9) + stderror + robust + endCell + beginCell + fmt(coe.getZ(), 5, 9) + endCell + beginCell + fmt(coe.getPvalue(), 6, 9) + endCell + beginCell + fmt(coe.getHazardRatio(), 3, 9) + endCell + beginCell + fmt(coe.getHazardRatioLoCI(), 3, 9) + endCell + beginCell + fmt(coe.getHazardRatioHiCI(), 3, 9) + endCell + endLine; 276 } 277 return o; 278 } 279 280 /** 281 * 282 * @param d 283 * @param precision 284 * @param pad 285 * @return 286 */ 287 public static String fmt(Double d, int precision, int pad) { 288 if(d == null) 289 return ""; 290 if(Double.isNaN(d)) 291 return ""; 292 String value = ""; 293 DecimalFormat dfe = new DecimalFormat("0.00E0"); 294 String dpad = "0."; 295 double p = 1.0; 296 for (int i = 0; i < (precision); i++) { 297 dpad = dpad + "0"; 298 p = p / 10.0; 299 } 300 DecimalFormat df = new DecimalFormat(dpad); 301 if (Math.abs(d) >= p) { 302 value = df.format(d); 303 } else { 304 value = dfe.format(d); 305 } 306 int length = value.length(); 307 int extra = pad - length; 308 if (extra > 0) { 309 for (int i = 0; i < extra; i++) { 310 value = " " + value; 311 } 312 } 313 return value; 314 } 315 316 /** 317 * 318 */ 319 private void calcSummaryValues() { 320 321 //beta 322 323 ArrayList<String> variables = new ArrayList<String>(coefficientsList.keySet()); 324 for (int i = 0; i < variables.size(); i++) { 325 String variable = variables.get(i); 326 CoxCoefficient coe = coefficientsList.get(variable); 327 coe.setStdError(Math.sqrt(imat[i][i])); //values can be updated to reflect new error 328 if (naive_imat != null) { 329 coe.setRobustStdError(Math.sqrt(naive_imat[i][i])); 330 } 331 coe.setZ(coe.getCoeff() / coe.getStdError()); 332 coe.setPvalue(ChiSq.norm(Math.abs(coe.getCoeff() / coe.getStdError()))); 333 //z <- qnorm((1 + conf.int)/2, 0, 1) 334 double z = 1.959964; 335 coe.setHazardRatioLoCI(Math.exp(coe.getCoeff() - z * coe.getStdError())); 336 coe.setHazardRatioHiCI(Math.exp(coe.getCoeff() + z * coe.getStdError())); 337 } 338 339 logTest = -2 * (loglikInit - loglikFinal); 340 logTestpval = ChiSq.chiSq(logTest, (int) degreeFreedom); 341 342 scoreLogrankTestpvalue = ChiSq.chiSq(scoreLogrankTest, (int) degreeFreedom); 343 if (rscore != null) { 344 rscoreLogrankTestpvalue = ChiSq.chiSq(rscore, (int) degreeFreedom); 345 } 346 } 347 348 /** 349 * 350 */ 351 public void dump() { 352 353 //need an ordered list for comparing to R dumps 354 355// ArrayList<SurvivalInfo> orderedSurvivalInfoList = new ArrayList<SurvivalInfo>(survivalInfoList); 356// SurvivalInfoComparator sicSort = new SurvivalInfoComparator(); 357 // Collections.sort(orderedSurvivalInfoList,sicSort); 358 359 360 System.out.println(); 361 System.out.println("$coef"); 362 for (CoxCoefficient coe : coefficientsList.values()) { 363 System.out.print(coe.getCoeff() + " "); 364 } 365 System.out.println(); 366 System.out.println("$means"); 367 368 for (CoxCoefficient coe : coefficientsList.values()) { 369 System.out.print(coe.getMean() + " "); 370 } 371 System.out.println(); 372 System.out.println("$u"); 373 374 for (double d : u) { 375 System.out.print(d + " "); 376 } 377 378 System.out.println(); 379 System.out.println("$imat"); 380 for (int i = 0; i < imat.length; i++) { 381 for (int j = 0; j < imat[0].length; j++) { 382 System.out.print(imat[i][j] + " "); 383 } 384 System.out.println(); 385 } 386 387 if (this.naive_imat != null) { 388 System.out.println("$naive_imat"); 389 for (int i = 0; i < naive_imat.length; i++) { 390 for (int j = 0; j < naive_imat[0].length; j++) { 391 System.out.print(naive_imat[i][j] + " "); 392 } 393 System.out.println(); 394 } 395 } 396 397 System.out.println(); 398 System.out.println("$loglik"); 399 400 System.out.println(loglikInit + " " + loglikFinal); 401 402 System.out.println(); 403 System.out.println("$sctest"); 404 405 System.out.println(this.scoreLogrankTest); 406 407 System.out.println("$iter"); 408 System.out.println(this.iterations); 409 410 System.out.println("$flag"); 411 System.out.println(flag); 412 413 System.out.println(); 414// if (false) { 415// System.out.println("ID LP Score Residuals"); 416// for (SurvivalInfo si : orderedSurvivalInfoList) { 417// System.out.println(si.getOrder() + " " + si.getLinearPredictor() + " " + si.getScore() + " " + si.getResidual()); 418// 419// } 420// System.out.println(); 421// ArrayList<String> variables = new ArrayList<String>(coefficientsList.keySet()); 422// System.out.print("Sample"); 423// for (String v : variables) { 424// System.out.print(" " + v); 425// } 426// System.out.println("rr"); 427// for (SurvivalInfo si : orderedSurvivalInfoList) { 428// System.out.print(si.getOrder()); 429// for (String v : variables) { 430// System.out.print(" " + si.getResidualVariable(v)); 431// } 432// System.out.println(); 433// } 434// } 435 436 } 437 438 /** 439 * 440 * @param d 441 * @param pad 442 * @return 443 */ 444 public String fmtpr(String d, int pad) { 445 int length = d.length(); 446 int extra = pad - length; 447 if (extra < 0) { 448 extra = 0; 449 } 450 String v = d; 451 for (int i = 0; i < extra; i++) { 452 v = v + " "; 453 } 454 return v; 455 } 456 457 /** 458 * Pad left a string with spaces 459 * 460 * @param d 461 * @param pad 462 * @return 463 */ 464 public String fmtpl(String d, int pad) { 465 int length = d.length(); 466 int extra = pad - length; 467 if (extra < 0) { 468 extra = 0; 469 } 470 String v = d; 471 for (int i = 0; i < extra; i++) { 472 v = " " + v; 473 } 474 return v; 475 } 476 477 @Override 478 public String toString() { 479 return toString("", " ", "\r\n"); 480 } 481 482 /** 483 * 484 * @param beginLine 485 * @param del 486 * @param endLine 487 * @return 488 */ 489 public String toString(String beginLine, String del, String endLine) { 490 491 492 493 String o = beginLine + fmtpl("", 9) + fmtpl("Avg", 9) + fmtpl("SD", 9) + endLine; 494 for (CoxCoefficient coe : coefficientsList.values()) { 495 o = o + beginLine + fmtpr(coe.getName(), 9) + fmt(coe.getMean(), 4, 9) + fmt(coe.getStandardDeviation(), 4, 9) + endLine; 496 497 } 498 499 o = o + beginLine + endLine; 500 501 502 503 o = o + beginLine + "n= " + this.numSamples + ", number of events=" + this.numEvents + endLine; 504 o = o + getCoefficientText(true, beginLine, del, "", endLine); 505 506 o = o + beginLine + endLine; 507 508 if (baselineSurvivorFunction.size() > 0) { 509 o = o + beginLine + "Baseline Survivor Function (at predictor means)" + endLine; 510 for (Double time : baselineSurvivorFunction.keySet()) { 511 Double mean = baselineSurvivorFunction.get(time); 512 o = o + beginLine + fmt(time, 4, 10) + fmt(mean, 4, 10) + endLine; 513 } 514 } 515 516 o = o + beginLine + endLine; 517 o = o + beginLine + "Overall Model Fit" + endLine; 518 o = o + beginLine + "Iterations=" + iterations + endLine; 519 520 521 o = o + beginLine + "Likelihood ratio test = " + fmt(this.logTest, 2, 0) + " df=" + this.degreeFreedom + " p-value=" + fmt(this.logTestpval, 7, 0) + endLine; 522 523 o = o + beginLine + "Wald test = " + fmt(waldTestInfo.getTest(), 2, 0) + " df=" + waldTestInfo.getDf() + " p-value=" + fmt(waldTestInfo.getPvalue(), 7, 0) + endLine; 524 o = o + beginLine + "Score (logrank) test = " + fmt(scoreLogrankTest, 2, 0) + " df=" + ((int) (this.degreeFreedom)) + " p-value=" + fmt(this.scoreLogrankTestpvalue, 7, 0); 525 526 if (this.rscore != null) { 527 o = o + ", Robust = " + fmt(rscore, 2, 0) + " p-value=" + fmt(rscoreLogrankTestpvalue, 7, 0); 528 529 } 530 531 o = o + endLine; 532 533 // o = o + "Rank of solution flag=" + flag + "\r\n"; 534 // o = o + "Log lik Initial=" + loglikInit + "\r\n"; 535 // o = o + "Log lik Final=" + loglikFinal + "\r\n"; 536 o = o + beginLine + "Method=" + method.name() + endLine; 537 538 539 540 541 return o; 542 } 543 544 /** 545 * @return the scoreLogrankTest 546 */ 547 public double getChiSquare() { 548 return scoreLogrankTest; 549 } 550 551 /** 552 * @return the degreeFreedom 553 */ 554 public double getDegreeFreedom() { 555 return degreeFreedom; 556 } 557 558 /** 559 * @return the scoreLogrankTestpvalue 560 */ 561 public double getOverallModelFitPvalue() { 562 return scoreLogrankTestpvalue; 563 } 564 565 /** 566 * @return the rscore 567 */ 568 public Double getRscore() { 569 return rscore; 570 } 571 572 /** 573 * @param rscore the rscore to set 574 */ 575 public void setRscore(Double rscore) { 576 this.rscore = rscore; 577 if (rscore != null) { 578 rscoreLogrankTestpvalue = ChiSq.chiSq(rscore, (int) degreeFreedom); 579 } 580 } 581 582 /** 583 * @return the rscoreLogrankTestpvalue 584 */ 585 public Double getRscoreLogrankTestpvalue() { 586 return rscoreLogrankTestpvalue; 587 } 588 589 /** 590 * @param rscoreLogrankTestpvalue the rscoreLogrankTestpvalue to set 591 */ 592 public void setRscoreLogrankTestpvalue(Double rscoreLogrankTestpvalue) { 593 this.rscoreLogrankTestpvalue = rscoreLogrankTestpvalue; 594 } 595 596 /** 597 * @return the scoreLogrankTest 598 */ 599 public Double getScoreLogrankTest() { 600 return scoreLogrankTest; 601 } 602 603 /** 604 * @param scoreLogrankTest the scoreLogrankTest to set 605 */ 606 public void setScoreLogrankTest(Double scoreLogrankTest) { 607 this.scoreLogrankTest = scoreLogrankTest; 608 } 609 610 /** 611 * @return the scoreLogrankTestpvalue 612 */ 613 public Double getScoreLogrankTestpvalue() { 614 return scoreLogrankTestpvalue; 615 } 616 617 /** 618 * @param scoreLogrankTestpvalue the scoreLogrankTestpvalue to set 619 */ 620 public void setScoreLogrankTestpvalue(Double scoreLogrankTestpvalue) { 621 this.scoreLogrankTestpvalue = scoreLogrankTestpvalue; 622 } 623 624 /** 625 * @return the metaDataFilter 626 */ 627 public LinkedHashMap<String, String> getMetaDataFilter() { 628 return metaDataFilter; 629 } 630 631 /** 632 * @param metaDataFilter the metaDataFilter to set 633 */ 634 public void setMetaDataFilter(LinkedHashMap<String, String> metaDataFilter) { 635 this.metaDataFilter = metaDataFilter; 636 } 637 638 /** 639 * @return the coefficientsList 640 */ 641 public LinkedHashMap<String, CoxCoefficient> getCoefficientsList() { 642 return coefficientsList; 643 } 644 645 /** 646 * @return the waldTestInfo 647 */ 648 public WaldTestInfo getWaldTestInfo() { 649 return waldTestInfo; 650 } 651 652 /** 653 * @return the imat 654 */ 655 public double[][] getImat() { 656 return imat; 657 } 658 659 /** 660 * @return the naive_imat 661 */ 662 public double[][] getNaive_imat() { 663 return naive_imat; 664 } 665 666 /** 667 * @param degreeFreedom the degreeFreedom to set 668 */ 669 public void setDegreeFreedom(double degreeFreedom) { 670 this.degreeFreedom = degreeFreedom; 671 } 672 673 /** 674 * @param waldTestInfo the waldTestInfo to set 675 */ 676 public void setWaldTestInfo(WaldTestInfo waldTestInfo) { 677 this.waldTestInfo = waldTestInfo; 678 } 679}