[Genabel-commits] r982 - in pkg/MixABEL: . R man src/MXlib
noreply at r-forge.r-project.org
noreply at r-forge.r-project.org
Tue Oct 23 19:12:51 CEST 2012
Author: yurii
Date: 2012-10-23 19:12:51 +0200 (Tue, 23 Oct 2012)
New Revision: 982
Added:
pkg/MixABEL/man/MixABEL-dash-package.Rd
Modified:
pkg/MixABEL/CHANGES.LOG
pkg/MixABEL/DESCRIPTION
pkg/MixABEL/R/fgls.R
pkg/MixABEL/R/zzz.R
pkg/MixABEL/man/FGLS.Rd
pkg/MixABEL/man/FastMixedModel.Rd
pkg/MixABEL/man/GWFGLS.Rd
pkg/MixABEL/man/MixABEL-package.Rd
pkg/MixABEL/man/summaryFGLS.Rd
pkg/MixABEL/src/MXlib/fglsWrapper.cpp
Log:
small changes + bug fix preventing gtcoding="probability" working (fix in line 523 of fgls.R)
Modified: pkg/MixABEL/CHANGES.LOG
===================================================================
--- pkg/MixABEL/CHANGES.LOG 2012-10-23 16:52:51 UTC (rev 981)
+++ pkg/MixABEL/CHANGES.LOG 2012-10-23 17:12:51 UTC (rev 982)
@@ -1,3 +1,7 @@
+***** v. 0.1-2 (2012.10.23)
+
+Fixed bug preventing correct results with gtcoding="probability"
+
***** v. 0.1-1 (2011.02.23)
Memory leak problem fix by William Astle
Modified: pkg/MixABEL/DESCRIPTION
===================================================================
--- pkg/MixABEL/DESCRIPTION 2012-10-23 16:52:51 UTC (rev 981)
+++ pkg/MixABEL/DESCRIPTION 2012-10-23 17:12:51 UTC (rev 982)
@@ -1,8 +1,8 @@
Package: MixABEL
Type: Package
Title: mixed models for genetic association analysis
-Version: 0.1-1
-Date: 2011-02-23
+Version: 0.1-2
+Date: 2012-10-23
Author: Yurii Aulchenko, William Astle, Erik Roos, Marcel Kempenaar
Maintainer: Yurii Aulchenko <yurii.aulchenko at gmail.com>
Depends: R (>= 2.4.0), MASS, mvtnorm, GenABEL (>= 1.5-8), DatABEL (>= 0.1-5)
Modified: pkg/MixABEL/R/fgls.R
===================================================================
--- pkg/MixABEL/R/fgls.R 2012-10-23 16:52:51 UTC (rev 981)
+++ pkg/MixABEL/R/fgls.R 2012-10-23 17:12:51 UTC (rev 982)
@@ -174,7 +174,7 @@
}
#' Genome-wide FGLS
#'
-#' Genome-wide Feasible GLS
+#' Genome-wide Feasible GLS
#'
#' @param formula analysis formula; should contain special 'SNP' term
#' @param data phenotypic data frame, or 'gwaa.data-class' object
@@ -384,7 +384,7 @@
}
rm(Xm,YmXB,bt,Y);gc()
}
- #print(sigma2.score)
+ #cat("sigma2.score =",sigma2.score,"\n")
if (!is.null(W)) {
require(MASS)
@@ -520,7 +520,7 @@
as.integer(gtNrow), as.integer(gtNcol),
as.character("fgls"),
"R", as.integer(2),
- as.integer(1), # STEP
+ as.integer(step), # STEP
as.integer(10),
SEXP_ptr_to_gsl_Y, SEXP_ptr_to_gsl_X,
SEXP_ptr_to_gsl_tXW_fixed, SEXP_ptr_to_gsl_W,
Modified: pkg/MixABEL/R/zzz.R
===================================================================
--- pkg/MixABEL/R/zzz.R 2012-10-23 16:52:51 UTC (rev 981)
+++ pkg/MixABEL/R/zzz.R 2012-10-23 17:12:51 UTC (rev 982)
@@ -1,3 +1,3 @@
.onLoad <- function(lib, pkg) {
- cat("MixABEL v 0.1-1 (February 23, 2011) loaded\n")
+ cat("MixABEL v 0.1-2 (Novenber 23, 2012) loaded\n")
}
\ No newline at end of file
Modified: pkg/MixABEL/man/FGLS.Rd
===================================================================
--- pkg/MixABEL/man/FGLS.Rd 2012-10-23 16:52:51 UTC (rev 981)
+++ pkg/MixABEL/man/FGLS.Rd 2012-10-23 17:12:51 UTC (rev 982)
@@ -1,21 +1,41 @@
\name{FGLS}
\alias{FGLS}
-\title{Feasible GLS...}
-\usage{FGLS(Y, X, W, test="wald", whichtest=c(FALSE, rep(TRUE, dim(X)[2] -
- 1)))}
-\description{Feasible GLS}
-\details{Feasible Generalised Least Squares}
-\value{List with elements 'beta' -- estimates fo the regression
-coefficients; 'V' -- variance covariance matrix for parameters
-estimates; 'T2' -- test statistics (distributed as Chi-squared under
-the null) for the testing of whichtest parameters; 'df' -- the number of
-degrees of freedom of the T2 test; 'tested' -- which parameters were
-tested with T2; 'meanColX' -- mean value of variable in columns of X;
-'n' -- length of Y (== height of X)}
-\author{Yurii Aulchenko}
-\arguments{\item{Y}{dependent variable}
-\item{X}{design matrix (including intercept, if necessary)}
-\item{test}{test to be applied, one of 'wald', 'score' or 'robust'}
-\item{whichtest}{which independent variables to be tested (set to 'TRUE')}
-\item{W}{for GLS, inverse variance-covariance matrix, as such returned by
-GenABEL's polygenic(...)$InvSigma, or NULL for LS}}
+\title{Feasible GLS}
+\usage{
+ FGLS(Y, X, W = NULL, test = "wald",
+ whichtest = c(FALSE, rep(TRUE, dim(X)[2] - 1)))
+}
+\arguments{
+ \item{Y}{dependent variable}
+
+ \item{X}{design matrix (including intercept, if
+ necessary)}
+
+ \item{test}{test to be applied, one of 'wald', 'score' or
+ 'robust'}
+
+ \item{whichtest}{which independent variables to be tested
+ (set to 'TRUE')}
+
+ \item{W}{for GLS, inverse variance-covariance matrix, as
+ such returned by GenABEL's polygenic(...)$InvSigma, or
+ NULL for LS}
+}
+\value{
+ List with elements 'beta' -- estimates fo the regression
+ coefficients; 'V' -- variance covariance matrix for
+ parameters estimates; 'T2' -- test statistics
+ (distributed as Chi-squared under the null) for the
+ testing of whichtest parameters; 'df' -- the number of
+ degrees of freedom of the T2 test; 'tested' -- which
+ parameters were tested with T2; 'meanColX' -- mean value
+ of variable in columns of X; 'n' -- length of Y (==
+ height of X)
+}
+\description{
+ Feasible Generalised Least Squares
+}
+\author{
+ Yurii Aulchenko
+}
+
Modified: pkg/MixABEL/man/FastMixedModel.Rd
===================================================================
--- pkg/MixABEL/man/FastMixedModel.Rd 2012-10-23 16:52:51 UTC (rev 981)
+++ pkg/MixABEL/man/FastMixedModel.Rd 2012-10-23 17:12:51 UTC (rev 982)
@@ -1,25 +1,38 @@
\name{FastMixedModel}
\alias{FastMixedModel}
-\title{fast mixed models...}
-\usage{FastMixedModel(Response, Explan, Kin, Covariates, nu_naught=0,
- gamma_naught=0)}
-\description{fast mixed models}
-\details{fast mixed models -- BETA VERSION
-If compiled against OMP this library can exploit multi-core parallelism
-Does not cope with missing data at present}
-\value{a list with values ...}
-\seealso{mmscore}
-\references{reference to fill in}
-\author{William Astle \email{fio at where}}
-\arguments{\item{Response}{is an n dimensional vector of Responses}
-\item{Explan}{is an n*p matrix of Explanatory variables, each to be tested marginally (SNPS)}
-\item{Kin}{is the n*n kinship matrix}
-\item{Covariates}{is an n*k matrix of Covariates}
-\item{nu_naught}{(and gamma_naught) are hyperparameters which control the heaviness
-of the tails of the test distribution (recommend leave them unchanged).}
-\item{gamma_naught}{(and nu_naught) are hyperparameters which control the heaviness
-of the tails of the test distribution (recommend leave them unchanged).}}
-\examples{require(mvtnorm)
+\title{fast mixed models}
+\usage{
+ FastMixedModel(Response, Explan, Kin, Covariates = NULL,
+ nu_naught = 0, gamma_naught = 0)
+}
+\arguments{
+ \item{Response}{is an n dimensional vector of Responses}
+
+ \item{Explan}{is an n*p matrix of Explanatory variables,
+ each to be tested marginally (SNPS)}
+
+ \item{Kin}{is the n*n kinship matrix}
+
+ \item{Covariates}{is an n*k matrix of Covariates}
+
+ \item{nu_naught}{(and gamma_naught) are hyperparameters
+ which control the heaviness of the tails of the test
+ distribution (recommend leave them unchanged).}
+
+ \item{gamma_naught}{(and nu_naught) are hyperparameters
+ which control the heaviness of the tails of the test
+ distribution (recommend leave them unchanged).}
+}
+\value{
+ a list with values ...
+}
+\description{
+ fast mixed models -- BETA VERSION If compiled against OMP
+ this library can exploit multi-core parallelism Does not
+ cope with missing data at present
+}
+\examples{
+require(mvtnorm)
data(ge03d2.clean)
df <- ge03d2.clean[1:250,autosomal(ge03d2.clean)]
NSNPS <- nsnps(df)
@@ -27,7 +40,7 @@
gkin <- ibs(df[,autosomal(df)],w="freq")
ngkin <- gkin
-ngkin[upper.tri(ngkin)] <- t(ngkin)[upper.tri(ngkin)]
+ngkin[upper.tri(ngkin)] <- t(ngkin)[upper.tri(ngkin)]
ngkin[1:5,1:5]
mysig <- (modh2*2*ngkin+(1.-modh2)*diag(dim(ngkin)[1]))
mysig[1:5,1:5]
@@ -53,13 +66,13 @@
expl <- as.numeric(df[,1:NSNPS])
summary(res)
covariates <- matrix(c(phdata(df)$sex,phdata(df)$age),ncol=2)
-summary(covariates)
+summary(covariates)
time0.fmm <- proc.time()
fmm <- FastMixedModel(Response=res,
-Explan=expl,
-Kin = gkin,
-Cov=covariates)
+ Explan=expl,
+ Kin = gkin,
+ Cov=covariates)
time.fmm <- proc.time() - time0.fmm
time.h2
@@ -73,4 +86,15 @@
fmm$null.herit
cor(mms[,"chi2.1df"],fmm$chi.sq)^2
-plot(mms[,"chi2.1df"],fmm$chi.sq)}
+plot(mms[,"chi2.1df"],fmm$chi.sq)
+}
+\author{
+ William Astle \email{fio at where}
+}
+\references{
+ reference to fill in
+}
+\seealso{
+ mmscore
+}
+
Modified: pkg/MixABEL/man/GWFGLS.Rd
===================================================================
--- pkg/MixABEL/man/GWFGLS.Rd 2012-10-23 16:52:51 UTC (rev 981)
+++ pkg/MixABEL/man/GWFGLS.Rd 2012-10-23 17:12:51 UTC (rev 982)
@@ -1,49 +1,87 @@
\name{GWFGLS}
\alias{GWFGLS}
-\title{Genome-wide FGLS...}
-\usage{GWFGLS(formula, data, subset, weights, na.action, contrasts, offset, W,
- inverse=TRUE, na.SNP="impute", mincall=0.95, residuals=FALSE,
- test="wald", model.SNP="additive", genodata, gtcoding="typed",
- verbosity=1, varcov=FALSE, include.means=TRUE, singular="ignore",
- with.lm=FALSE, old=FALSE)}
-\description{Genome-wide FGLS}
-\details{Genome-wide Feasible GLS}
-\author{Yurii Aulchenko}
-\arguments{\item{formula}{analysis formula; should contain special 'SNP' term}
-\item{data}{phenotypic data frame, or 'gwaa.data-class' object}
-\item{subset}{subset of data (individuals)}
-\item{weights}{RESERVED FOR FURTURE USE}
-\item{na.action}{RESERVED FOR FURTURE USE}
-\item{contrasts}{RESERVED FOR FURTURE USE}
-\item{offset}{RESERVED FOR FURTURE USE}
-\item{W}{for GLS, (inverse of) variance-covariance matrix, as such returned by
-GenABEL's polygenic(...)\$InvSigma, or NULL for LS}
-\item{inverse}{wheather W is already inverted}
-\item{na.SNP}{how to deal with missing SNP data; 'impute' -- substitute
-with mean, 'drop' -- drop rows with missing genotypes}
-\item{mincall}{minimall call rate for a SNP (if less, the SNP is dropped)}
-\item{residuals}{use residuals for analysis? (faster, but less precise)}
-\item{test}{test to be applied, one of 'wald', 'score' or 'robust'}
-\item{model.SNP}{SNP model to apply, one of
-c("additive","dominantB","recessiveB","overdominant","genotypic")}
-\item{genodata}{genotypic data; can be missing if data is of 'gwaa.data-class'.
-Otherwise can be regular matrix or 'databel' matrix}
-\item{gtcoding}{one of c("typed","dose","probability")
-'typed' -- coded with NA, 0, 1, 2}
-\item{verbosity}{what to report; 0 -- only test stats,
-1 -- test stats and estimates concerning 'whichtest' parameters,
-2 -- test stats and estimates of all parameters}
-\item{varcov}{wheather var-cov matrix for estimated parameters
-is to be reported}
-\item{include.means}{weather mean values of predictors are to
-be reported}
-\item{singular}{what to do with linear dependencies in X (now only
-'ignore' implemented)}
-\item{with.lm}{whether LM should be run along (only test purposes; always
-falls back to 'old' R-only implementation)}
-\item{old}{if TRUE, old R-only code implementation is running (testing purposes,
-slow)}}
-\examples{library(MASS)
+\title{Genome-wide FGLS}
+\usage{
+ GWFGLS(formula, data, subset, weights, na.action,
+ contrasts = NULL, offset, W = NULL, inverse = TRUE,
+ na.SNP = "impute", mincall = 0.95, residuals = FALSE,
+ test = "wald", model.SNP = "additive", genodata,
+ gtcoding = "typed", verbosity = 1, varcov = FALSE,
+ include.means = TRUE, singular = "ignore",
+ with.lm = FALSE, old = FALSE)
+}
+\arguments{
+ \item{formula}{analysis formula; should contain special
+ 'SNP' term}
+
+ \item{data}{phenotypic data frame, or 'gwaa.data-class'
+ object}
+
+ \item{subset}{subset of data (individuals)}
+
+ \item{weights}{RESERVED FOR FURTURE USE}
+
+ \item{na.action}{RESERVED FOR FURTURE USE}
+
+ \item{contrasts}{RESERVED FOR FURTURE USE}
+
+ \item{offset}{RESERVED FOR FURTURE USE}
+
+ \item{W}{for GLS, (inverse of) variance-covariance
+ matrix, as such returned by GenABEL's
+ polygenic(...)\$InvSigma, or NULL for LS}
+
+ \item{inverse}{wheather W is already inverted}
+
+ \item{na.SNP}{how to deal with missing SNP data; 'impute'
+ -- substitute with mean, 'drop' -- drop rows with missing
+ genotypes}
+
+ \item{mincall}{minimall call rate for a SNP (if less, the
+ SNP is dropped)}
+
+ \item{residuals}{use residuals for analysis? (faster, but
+ less precise)}
+
+ \item{test}{test to be applied, one of 'wald', 'score' or
+ 'robust'}
+
+ \item{model.SNP}{SNP model to apply, one of
+ c("additive","dominantB","recessiveB","overdominant","genotypic")}
+
+ \item{genodata}{genotypic data; can be missing if data is
+ of 'gwaa.data-class'. Otherwise can be regular matrix or
+ 'databel' matrix}
+
+ \item{gtcoding}{one of c("typed","dose","probability")
+ 'typed' -- coded with NA, 0, 1, 2}
+
+ \item{verbosity}{what to report; 0 -- only test stats, 1
+ -- test stats and estimates concerning 'whichtest'
+ parameters, 2 -- test stats and estimates of all
+ parameters}
+
+ \item{varcov}{wheather var-cov matrix for estimated
+ parameters is to be reported}
+
+ \item{include.means}{weather mean values of predictors
+ are to be reported}
+
+ \item{singular}{what to do with linear dependencies in X
+ (now only 'ignore' implemented)}
+
+ \item{with.lm}{whether LM should be run along (only test
+ purposes; always falls back to 'old' R-only
+ implementation)}
+
+ \item{old}{if TRUE, old R-only code implementation is
+ running (testing purposes, slow)}
+}
+\description{
+ Genome-wide Feasible GLS
+}
+\examples{
+library(MASS)
library(mvtnorm)
library(GenABEL)
data(ge03d2.clean)
@@ -63,7 +101,7 @@
grel <- gkin
grel[upper.tri(grel)] <- t(grel)[upper.tri(grel)]
grel <- 2*grel
-Y <- as.vector(rmvnorm(1,mean=(X \%*\% betas_cov),sigma=s2*(modelh2*grel+diag(dim(grel)[1])*(1-modelh2))))
+Y <- as.vector(rmvnorm(1,mean=(X \\\%*\\\% betas_cov),sigma=s2*(modelh2*grel+diag(dim(grel)[1])*(1-modelh2))))
length(Y)
Y[2] <- NA
df at phdata$Y <- Y
@@ -85,5 +123,10 @@
aIN
aWN
complex_model <- GWFGLS(Y~age+sex*SNP,data=phdata(df),W=h2$InvSigma,
-genodata = gtReal,verbosity=2,varcov=TRUE,include.means=TRUE,model="genotypic")
-complex_model}
+ genodata = gtReal,verbosity=2,varcov=TRUE,include.means=TRUE,model="genotypic")
+complex_model
+}
+\author{
+ Yurii Aulchenko
+}
+
Added: pkg/MixABEL/man/MixABEL-dash-package.Rd
===================================================================
--- pkg/MixABEL/man/MixABEL-dash-package.Rd (rev 0)
+++ pkg/MixABEL/man/MixABEL-dash-package.Rd 2012-10-23 17:12:51 UTC (rev 982)
@@ -0,0 +1,8 @@
+\name{MixABEL-package}
+\alias{MixABEL-package}
+\title{MixABEL package...}
+\description{MixABEL package}
+\details{Package for fast (genome-wide) association analysis
+using mixed models. Main functions:
+\link{FastMixedModel}, \link{GWFGLS}}
+\author{Yurii Aulchenko, William Astle, Erik Roos, Marcel Kempenaar}
Property changes on: pkg/MixABEL/man/MixABEL-dash-package.Rd
___________________________________________________________________
Added: svn:mime-type
+ text/plain
Modified: pkg/MixABEL/man/MixABEL-package.Rd
===================================================================
--- pkg/MixABEL/man/MixABEL-package.Rd 2012-10-23 16:52:51 UTC (rev 981)
+++ pkg/MixABEL/man/MixABEL-package.Rd 2012-10-23 17:12:51 UTC (rev 982)
@@ -1,8 +1,16 @@
\name{MixABEL-package}
\alias{MixABEL-package}
-\title{MixABEL package...}
-\description{MixABEL package}
-\details{Package for fast (genome-wide) association analysis
-using mixed models. Main functions:
-\link{FastMixedModel}, \link{GWFGLS}}
-\author{Yurii Aulchenko, William Astle, Erik Roos, Marcel Kempenaar}
+\title{MixABEL package}
+\usage{
+ MixABEL-package()
+}
+\description{
+ Package for fast (genome-wide) association analysis using
+ mixed models. Main functions: \link{FastMixedModel},
+ \link{GWFGLS}
+}
+\author{
+ Yurii Aulchenko, William Astle, Erik Roos, Marcel
+ Kempenaar
+}
+
Modified: pkg/MixABEL/man/summaryFGLS.Rd
===================================================================
--- pkg/MixABEL/man/summaryFGLS.Rd 2012-10-23 16:52:51 UTC (rev 981)
+++ pkg/MixABEL/man/summaryFGLS.Rd 2012-10-23 17:12:51 UTC (rev 982)
@@ -1,16 +1,31 @@
\name{summaryFGLS}
\alias{summaryFGLS}
-\title{summary for FGLS...}
-\usage{summaryFGLS(x, verbosity=1, varcov=FALSE, include.means=FALSE)}
-\description{summary for FGLS}
-\details{summary for Feasible GLS}
-\value{one-line summary of output from FGLS function}
-\author{Yurii Aulchenko}
-\arguments{\item{x}{object returned by FGLS function}
-\item{verbosity}{what to report; 0 -- only test stats,
-1 -- test stats and estimates concerning 'whichtest' parameters,
-2 -- test stats and estimates of all parameters}
-\item{varcov}{wheather var-cov matrix for estimated parameters
-is to be reported}
-\item{include.means}{weather mean values of predictors are to
-be reported}}
+\title{summary for FGLS}
+\usage{
+ summaryFGLS(x, verbosity = 1, varcov = FALSE,
+ include.means = FALSE)
+}
+\arguments{
+ \item{x}{object returned by FGLS function}
+
+ \item{verbosity}{what to report; 0 -- only test stats, 1
+ -- test stats and estimates concerning 'whichtest'
+ parameters, 2 -- test stats and estimates of all
+ parameters}
+
+ \item{varcov}{wheather var-cov matrix for estimated
+ parameters is to be reported}
+
+ \item{include.means}{weather mean values of predictors
+ are to be reported}
+}
+\value{
+ one-line summary of output from FGLS function
+}
+\description{
+ summary for Feasible GLS
+}
+\author{
+ Yurii Aulchenko
+}
+
Modified: pkg/MixABEL/src/MXlib/fglsWrapper.cpp
===================================================================
--- pkg/MixABEL/src/MXlib/fglsWrapper.cpp 2012-10-23 16:52:51 UTC (rev 981)
+++ pkg/MixABEL/src/MXlib/fglsWrapper.cpp 2012-10-23 17:12:51 UTC (rev 982)
@@ -262,7 +262,7 @@
**/
// copy iterated data
- for (unsigned int i = 0; i < ncol; i++)
+ for (unsigned int i = 0; i < ncol; i++) {
for (unsigned int j = 0; j < nrow; j++) {
if (INTEGER(argList[5])[i] == 0) {}
else if (INTEGER(argList[5])[i] == 1) {
@@ -278,8 +278,21 @@
return;
}
}
+ }
delete [] processedGT;
+/** printf("nGTcols = %d\n",nGTcols);
+ printf("Matrix m\n");
+ for (int i=0;i<nrow;i++)
+ {
+ for (int j=0;j<ncol;j++)
+ {
+ printf("%f ",gsl_matrix_get(gsl_X,i,j));
+ }
+ printf("\n");
+ }
+ printf("\n");
+**/
// copy iterated data
//for (unsigned int i = noldcol; i < ncol; i++)
// for (unsigned int j = 0; j < nrow; j++)
More information about the Genabel-commits
mailing list