[Genabel-commits] r979 - branches/ProbABEL-refactoring/ProbABEL/src

[noreply at r-forge.r-project.org]

Author: maartenk
Date: 2012-10-14 14:59:36 +0200 (Sun, 14 Oct 2012)
New Revision: 979

Modified:
   branches/ProbABEL-refactoring/ProbABEL/src/maskedmatrix.cpp
   branches/ProbABEL-refactoring/ProbABEL/src/reg1.h
Log:
revert svn commit 852  (now the invarmatrix is passed by reference instead of value)
fixed bug in maskedmarix: after initialization the mask_of_old pointed to newmask and the masked matrix was not updated

Modified: branches/ProbABEL-refactoring/ProbABEL/src/maskedmatrix.cpp
===================================================================
--- branches/ProbABEL-refactoring/ProbABEL/src/maskedmatrix.cpp	2012-10-11 14:32:29 UTC (rev 978)
+++ branches/ProbABEL-refactoring/ProbABEL/src/maskedmatrix.cpp	2012-10-14 12:59:36 UTC (rev 979)
@@ -51,7 +51,8 @@
 masked_matrix::~masked_matrix() {
 	// TODO Auto-generated destructor stub
 }
-void masked_matrix::update_mask( short unsigned int *newmask) {
+void masked_matrix::update_mask(short unsigned int *newmask) {
+
 //find length of masked matrix
 	int nmeasured = 0;
 	for (int i = 0; i < length_of_mask; i++) {
@@ -70,7 +71,11 @@
 			masked_data = &matrix_masked_data;
 		} else {
 			//new mask differs from old matrix and create new.
-			mask_of_old = newmask;
+//			//mask_of_old = newmask;
+			//TODO(maarten): there must be a smarter way to copy these values
+			for (int i = 0; i < length_of_mask; i++) {
+				mask_of_old[i] = newmask[i];
+			}
 			mask_symmetric(nmeasured);
 			masked_data = &matrix_masked_data;
 		}
@@ -93,8 +98,6 @@
 				}
 			i1++;
 		}
-
-
 }
 //mematrix<double>* masked_matrix::get_matrix()
 //{

Modified: branches/ProbABEL-refactoring/ProbABEL/src/reg1.h
===================================================================
--- branches/ProbABEL-refactoring/ProbABEL/src/reg1.h	2012-10-11 14:32:29 UTC (rev 978)
+++ branches/ProbABEL-refactoring/ProbABEL/src/reg1.h	2012-10-14 12:59:36 UTC (rev 979)
@@ -243,7 +243,7 @@
 
 	void base_score(mematrix<double> &resid, regdata &rdata, int verbose,
 			double tol_chol, int model, int interaction, int ngpreds,
-			masked_matrix invvarmatrix, int nullmodel = 0) {
+			const masked_matrix &invvarmatrix, int nullmodel = 0) {
 		mematrix<double> oX = rdata.extract_genotypes();
 		mematrix<double> X = apply_model(oX, model, interaction, ngpreds,
 				rdata.is_interaction_excluded, false, nullmodel);
@@ -334,7 +334,7 @@
 //    }
 
 	void estimate(regdata& rdatain, int verbose, double tol_chol, int model,
-			int interaction, int ngpreds, masked_matrix invvarmatrixin,
+			int interaction, int ngpreds, masked_matrix &invvarmatrixin,
 			int robust, int nullmodel = 0) {
 		//suda ineraction parameter
 		// model should come here
@@ -568,7 +568,7 @@
 
 	void score(mematrix<double> &resid, regdata &rdatain, int verbose,
 			double tol_chol, int model, int interaction, int ngpreds,
-			masked_matrix invvarmatrix, int nullmodel = 0) {
+			const masked_matrix &invvarmatrix, int nullmodel = 0) {
 		regdata rdata = rdatain.get_unmasked_data();
 		base_score(resid, rdata, verbose, tol_chol, model, interaction, ngpreds,
 				invvarmatrix, nullmodel = 0);
@@ -804,64 +804,3 @@
 				invvarmatrix, nullmodel = 0);
 	}
 };
-
-//class coxph_reg
-//{
-//public:
-//    mematrix<double> beta;
-//    mematrix<double> sebeta;
-//    mematrix<double> residuals;
-//    double sigma2;
-//    double loglik;
-//    double chi2_score;
-//    int niter;
-//
-//    coxph_reg(coxph_data &cdatain)
-//    {
-//        coxph_data cdata = cdatain.get_unmasked_data();
-//        beta.reinit(cdata.X.nrow, 1);
-//        sebeta.reinit(cdata.X.nrow, 1);
-//        loglik = -9.999e+32;
-//        sigma2 = -1.;
-//        chi2_score = -1.;
-//        niter = 0;
-//    }
-//    ~coxph_reg()
-//    {
-//        //		delete beta;
-//        //		delete sebeta;
-//    }
-//    void estimate(coxph_data &cdatain, int verbose, int maxiter, double eps,
-//            double tol_chol, int model, int interaction, int ngpreds,
-//            bool iscox, int nullmodel = 0)
-//    {
-//        //		cout << "model = " << model << "\n";
-//        //		cdata.X.print();
-//        coxph_data cdata = cdatain.get_unmasked_data();
-//        mematrix<double> X = t_apply_model(cdata.X, model, interaction, ngpreds,
-//                iscox, nullmodel);
-//        //		X.print();
-//        int length_beta = X.nrow;
-//        beta.reinit(length_beta, 1);
-//        sebeta.reinit(length_beta, 1);
-//        mematrix<double> newoffset = cdata.offset;
-//        newoffset = cdata.offset - (cdata.offset).column_mean(0);
-//        mematrix<double> means(X.nrow, 1);
-//        for (int i = 0; i < X.nrow; i++)
-//            beta[i] = 0.;
-//        mematrix<double> u(X.nrow, 1);
-//        mematrix<double> imat(X.nrow, X.nrow);
-//        double work[X.ncol * 2 + 2 * (X.nrow) * (X.nrow) + 3 * (X.nrow)];
-//        double loglik_int[2];
-//        int flag;
-//        double sctest = 1.0;
-////TODO(maarten): remove the following comment signs. This is done for testing purpose only EVIL!
-//        coxfit2(&maxiter, &cdata.nids, &X.nrow, cdata.stime.data,
-//                cdata.sstat.data, X.data, newoffset.data, cdata.weights.data,
-//                cdata.strata.data, means.data, beta.data, u.data, imat.data,loglik_int, &flag, work, &eps, &tol_chol, &sctest);
-//        for (int i = 0; i < X.nrow; i++)
-//            sebeta[i] = sqrt(imat.get(i, i));
-//        loglik = loglik_int[1];
-//        niter = maxiter;
-//    }
-//};