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.Cholesky2;
024import org.biojava.nbio.survival.cox.stats.Chsolve2;
025
026/**
027 *
028 * @author Scooter Willis <willishf at gmail dot com>
029 */
030public class WaldTest {
031//coxph_wtest, df=as.integer(nvar),as.integer(ntest),as.double(var),tests= as.double(b),solve= double(nvar*ntest),as.double(toler.chol))
032        //coxph_wtest(Sint *nvar2, Sint *ntest, double *var, double *b,double *solve, double *tolerch)
033
034        /**
035         *
036         * @param var
037         * @param b
038         * @param toler_chol
039         * @return
040         */
041        public static WaldTestInfo process(double[][] var, double[] b, double toler_chol) {
042                double[][] b_ = new double[1][b.length];
043
044                for(int i = 0; i < b.length; i++){
045                        b_[0][i] = b[i];
046                }
047
048                return process(var,b_,toler_chol);
049
050        }
051
052        /**
053         *
054         * @param var
055         * @param b
056         * @param toler_chol
057         * @return
058         */
059        public static WaldTestInfo process(double[][] var, double[][] b, double toler_chol) {
060
061
062                int i = 0;
063
064        //      if(ci.coefficientsList.size() == 1){
065        //          double b_ = b[0][i];
066        //          double t = (b_ * b_) / var[0][0];
067        //          return;
068        //      }
069
070
071                //  double toler_chol = ci.toler;
072                int ntest = 1;
073                int nvar = b[0].length;
074                double sum = 0;
075                double[][] solve = new double[ntest][nvar];
076                double[] bsum = new double[ntest];
077
078                Cholesky2.process(var, nvar, toler_chol);
079
080                int df = 0;
081                for (i = 0; i < nvar; i++) {
082                        if (var[i][i] > 0) {
083                                df++;  /* count up the df */
084                        }
085                }
086
087                for (i = 0; i < ntest; i++) {
088                        for (int j = 0; j < nvar; j++) {
089                                solve[i][j] = b[i][j];
090                        }
091                        Chsolve2.process(var, nvar, solve, i);   /*solve now has b* var-inverse */
092
093                        sum = 0;
094                        for (int j = 0; j < nvar; j++) {
095                                sum += b[i][j] * solve[i][j];
096                        }
097                        bsum[i] = sum;                     /* save the result */
098                        //b += nvar;    /*move to next column of b */
099                        // solve += nvar;
100                }
101                //* nvar2 = df;
102                WaldTestInfo waldTestInfo = new WaldTestInfo();
103
104                waldTestInfo.setDf(df);
105                waldTestInfo.solve = solve;
106                waldTestInfo.bsum = bsum;
107
108                return waldTestInfo;
109        }
110}
111