[Genabel-commits] r1104 - in branches/ProbABEL-pacox/v.0.3.0/ProbABEL: doc src
noreply at r-forge.r-project.org
noreply at r-forge.r-project.org
Wed Feb 13 22:06:37 CET 2013
Author: lckarssen
Date: 2013-02-13 22:06:37 +0100 (Wed, 13 Feb 2013)
New Revision: 1104
Modified:
branches/ProbABEL-pacox/v.0.3.0/ProbABEL/doc/pacoxph.1
branches/ProbABEL-pacox/v.0.3.0/ProbABEL/src/gendata.h
branches/ProbABEL-pacox/v.0.3.0/ProbABEL/src/main.cpp
branches/ProbABEL-pacox/v.0.3.0/ProbABEL/src/reg1.cpp
Log:
- reg1.cpp: IMPORTANT STUFF: commented the interaction_only if statement. This makes the probability input work, but obviously breaks functionality. THIS NEEDS TO BE FIXED before merging back to trunk.
- main.cpp: Minor code layout improvement, but also still contains some commented debug prints
- gendata.h: removed commented statement
- man page: remved message about pacoxph being buggy
Modified: branches/ProbABEL-pacox/v.0.3.0/ProbABEL/doc/pacoxph.1
===================================================================
--- branches/ProbABEL-pacox/v.0.3.0/ProbABEL/doc/pacoxph.1 2013-02-13 09:36:42 UTC (rev 1103)
+++ branches/ProbABEL-pacox/v.0.3.0/ProbABEL/doc/pacoxph.1 2013-02-13 21:06:37 UTC (rev 1104)
@@ -1,4 +1,4 @@
-.TH pacoxph 1 "23 February 2012"
+.TH pacoxph 1 "08 February 2013"
.SH NAME
pacoxph \- Perform Genome-Wide Association Analysis using a linear model
.SH SYNOPSIS
@@ -78,8 +78,6 @@
.SH "SEE ALSO"
palinear(1), palogist(1)
.SH BUGS
-Unfortunately
-.B pacoxph
-is in a buggy state at the moment. It cannot use files in DatABEL format.
+None known at the moment. The bugtracker is located at https://r-forge.r-project.org/tracker/?group_id=505
.SH AUTHORS
Lennart C. Karssen
Modified: branches/ProbABEL-pacox/v.0.3.0/ProbABEL/src/gendata.h
===================================================================
--- branches/ProbABEL-pacox/v.0.3.0/ProbABEL/src/gendata.h 2013-02-13 09:36:42 UTC (rev 1103)
+++ branches/ProbABEL-pacox/v.0.3.0/ProbABEL/src/gendata.h 2013-02-13 21:06:37 UTC (rev 1104)
@@ -43,7 +43,6 @@
mematrix<float> G;
AbstractMatrix * DAG;
unsigned short int * DAGmask;
- // mematrix<double> G;
};
#endif /* GENDATA_H_ */
Modified: branches/ProbABEL-pacox/v.0.3.0/ProbABEL/src/main.cpp
===================================================================
--- branches/ProbABEL-pacox/v.0.3.0/ProbABEL/src/main.cpp 2013-02-13 09:36:42 UTC (rev 1103)
+++ branches/ProbABEL-pacox/v.0.3.0/ProbABEL/src/main.cpp 2013-02-13 21:06:37 UTC (rev 1104)
@@ -359,8 +359,8 @@
nrd.estimate(nrgd, 0, CHOLTOL, 0, input_var.getInteraction(),
input_var.getNgpreds(), invvarmatrix, input_var.getRobust(), 1);
#elif COXPH
- coxph_reg nrd(nrgd);
-
+ coxph_reg nrd = coxph_reg(nrgd);
+// LCK std::cout << "Starting estimation of null model...\n";
nrd.estimate(nrgd, 0, MAXITER, EPS, CHOLTOL, 0,
input_var.getInteraction(), input_var.getNgpreds(), 1);
#endif
@@ -372,8 +372,8 @@
#else
regdata rgd(phd, gtd, 0, input_var.isIsInteractionExcluded());
#endif
-
std::cout << " formed regression object ...";
+// rgd.X.print();
std::cout << " done\n" << std::flush;
//________________________________________________________________
@@ -498,11 +498,13 @@
// freq = (gtd.G).column_mean(csnp)/2.;
gtd.get_var(csnp, snpdata1);
for (unsigned int ii = 0; ii < gtd.nids; ii++)
+ {
if (!isnan(snpdata1[ii]))
{
gcount++;
freq += snpdata1[ii] * 0.5;
}
+ }
}
freq /= static_cast<double>(gcount);
int poly = 1;
@@ -535,7 +537,7 @@
for (int model = 0; model < maxmod; model++)
{
- if (poly) //allel freq is not to rare
+ if (poly) //allele freq is not too rare
{
#if LOGISTIC
logistic_reg rd(rgd);
@@ -567,9 +569,11 @@
}
#elif COXPH
coxph_reg rd(rgd);
+// LCK std::cout << "HERE, model="<<model<<"\n";
rd.estimate(rgd, 0, MAXITER, EPS, CHOLTOL, model,
input_var.getInteraction(), true,
input_var.getNgpreds());
+// LCK std::cout << "HERE DONE\n";
#endif
if (!input_var.getAllcov() && model == 0
Modified: branches/ProbABEL-pacox/v.0.3.0/ProbABEL/src/reg1.cpp
===================================================================
--- branches/ProbABEL-pacox/v.0.3.0/ProbABEL/src/reg1.cpp 2013-02-13 09:36:42 UTC (rev 1103)
+++ branches/ProbABEL-pacox/v.0.3.0/ProbABEL/src/reg1.cpp 2013-02-13 21:06:37 UTC (rev 1104)
@@ -150,60 +150,65 @@
nX.reinit(X.nrow, (X.ncol - 1));
}
- int c1 = X.ncol - 2;
- int c2 = X.ncol - 1;
+ int c1 = X.ncol - 2; // column with Prob(A1A2)
+ int c2 = X.ncol - 1; // column with Prob(A1A1). Note the order is swapped cf the file!
for (int i = 0; i < X.nrow; i++)
for (int j = 0; j < (X.ncol - 2); j++)
nX[i * nX.ncol + j] = X[i * X.ncol + j];
for (int i = 0; i < nX.nrow; i++)
{
- if (model == 1)
+ if (model == 1) // additive
nX[i * nX.ncol + c1] = X[i * X.ncol + c1] + 2. * X[i * X.ncol + c2];
- else if (model == 2)
+ else if (model == 2) //dominant
nX[i * nX.ncol + c1] = X[i * X.ncol + c1] + X[i * X.ncol + c2];
- else if (model == 3)
+ else if (model == 3) // recessive
nX[i * nX.ncol + c1] = X[i * X.ncol + c2];
- else if (model == 4)
+ else if (model == 4) // over-dominant
nX[i * nX.ncol + c1] = X[i * X.ncol + c1];
if (interaction != 0)
nX[i * nX.ncol + c2] = X[i * nX.ncol + interaction]
* nX[i * nX.ncol + c1]; //Maksim: interaction with SNP
}
+
//Han Chen
- if (is_interaction_excluded)
- {
- mematrix<double> nX_without_interact_phe;
- nX_without_interact_phe.reinit(nX.nrow, nX.ncol - 1);
- int col_new;
- for (int row = 0; row < nX.nrow; row++)
- {
- col_new = -1;
- for (int col = 0; col < nX.ncol; col++)
- {
- if (col != interaction && !iscox)
- {
- col_new++;
- nX_without_interact_phe[row * nX_without_interact_phe.ncol
- + col_new] = nX[row * nX.ncol + col];
- }
- if (col != interaction - 1 && iscox)
- {
- col_new++;
- nX_without_interact_phe[row * nX_without_interact_phe.ncol
- + col_new] = nX[row * nX.ncol + col];
- }
- }
- }
- return nX_without_interact_phe;
- } //interaction_only, model!=0, ngpreds==2
+ // if (is_interaction_excluded)
+ // {
+ // mematrix<double> nX_without_interact_phe;
+ // nX_without_interact_phe.reinit(nX.nrow, nX.ncol - 1);
+ // int col_new;
+ // for (int row = 0; row < nX.nrow; row++)
+ // {
+ // col_new = -1;
+ // for (int col = 0; col < nX.ncol; col++)
+ // {
+ // if (col != interaction && !iscox)
+ // {
+ // col_new++;
+ // nX_without_interact_phe[row * nX_without_interact_phe.ncol
+ // + col_new] = nX[row * nX.ncol + col];
+ // }
+ // if (col != interaction - 1 && iscox)
+ // {
+ // col_new++;
+ // nX_without_interact_phe[row * nX_without_interact_phe.ncol
+ // + col_new] = nX[row * nX.ncol + col];
+ // }
+ // }
+ // }
+ // std::cout<<"LEAVINGHERE\n";
+ // return nX_without_interact_phe;
+ // } //interaction_only, model!=0, ngpreds==2
+// LCK printf("In apply_model: nX.nrow: %i nX.ncol: %i\n", nX.nrow, nX.ncol);
return nX;
}
mematrix<double> t_apply_model(mematrix<double>& X, int model, int interaction,
int ngpreds, bool iscox, int nullmodel)
{
+// LCK std::cout << "t_apply: ngpreds: " << ngpreds << "; model: "<< model<<"\n";
+
mematrix<double> tmpX = transpose(X);
mematrix<double> nX = apply_model(tmpX, model, interaction, ngpreds, iscox,
nullmodel);
More information about the Genabel-commits
mailing list