Hi,
I know the data frame c.df has the variables named exactly the way I want it
to be. I tried reading each covariate seperately but still it gives me an
error. I tried fidning out the error but it doesnt seem very evident. Here is
the error message
> #######################
>
>
> #retreive data
> # considering four covariates
> d.df=as.data.frame(read.table("c:/tina/phd/thesis/data/modified_data1.1.txt",header=T,sep=","))
> y=d.df[,ncol(d.df)]
> x=d.df[,1:4]
> # read each column of x seperately
> x1=d.df[,1]
> x2=d.df[,2]
> x3=d.df[,3]
> x4=d.df[,4]
> c.df=cbind(y,x)
> print(c.df)
y X1 X2 X3 X4
1 1 2 0 0 1
2 1 10 0 0 1
3 0 4 0 0 1
4 1 0 0 1 1
5 1 0 0 1 1
6 1 7 0 0 1
7 1 1 0 0 1
8 0 0 0 0 1
9 1 0 0 0 1
10 1 5 0 0 1
> p <- ncol(c.df)
>
> # setting error handler to true
> options(show.error.mesages = TRUE)
>
> # marginal log-prior of beta[]
> logpriorfun <- function(beta, mu, gshape, grate)
+ {
+ logprior = -p*log(2) + log(gamma(p+gshape)) - log(gamma(gshape))
+ + gshape*log(grate) - (p+gshape)* log(grate+sum(abs(beta)))
+ return(logprior)
+ }
> require(MCMCpack)
Loading required package: MCMCpack
Loading required package: coda
Loading required package: lattice
Loading required package: MASS
##
## Markov Chain Monte Carlo Package (MCMCpack)
## Copyright (C) 2003-2007 Andrew D. Martin and Kevin M. Quinn
##
## Support provided by the U.S. National Science Foundation
## (Grants SES-0350646 and SES-0350613)
##
[1] TRUE
Warning message:
package 'MASS' was built under R version 2.4.1
> a0 = 0.5
> b0 = 1
> mu0 = 0
> beta.init=list(c(0, rep(0.1,4)), c(0, rep(-0.1,4)), c(0, rep(0, 4)))
> burnin.cycles = 1000
> mcmc.cycles = 25000
> # three chains
> post.list <- lapply(beta.init, function(vec)
+ {
+ posterior <- MCMClogit(y~x1+x2+x3+x4, data=c.df, burnin=burnin.cycles,
mcmc=mcmc.cycles,
+ thin=5, tune=0.5, beta.start=vec, user.prior.density=logpriorfun, logfun=TRUE,
+ mu=mu0, gshape=a0, grate=b0)
+ return(posterior)
+ })
Error in tune %*% V : non-conformable arguments
>
> # for tracing error mesages
> geterrmessage()
[1] "Error in tune %*% V : non-conformable arguments\n"
> traceback()
3: MCMClogit(y ~ x1 + x2 + x3 + x4, data = c.df, burnin = burnin.cycles,
mcmc = mcmc.cycles, thin = 5, tune = 0.5, beta.start = vec,
user.prior.density = logpriorfun, logfun = TRUE, mu = mu0,
gshape = a0, grate = b0)
2: FUN(X[[1]], ...)
1: lapply(beta.init, function(vec) {
posterior <- MCMClogit(y ~ x1 + x2 + x3 + x4, data = c.df,
burnin = burnin.cycles, mcmc = mcmc.cycles, thin = 5,
tune = 0.5, beta.start = vec, user.prior.density = logpriorfun,
logfun = TRUE, mu = mu0, gshape = a0, grate = b0)
return(posterior)
})
>
Thanks for your help,
Anamika
Ravi Varadhan <[EMAIL PROTECTED]> wrote:
As the error message clearly indicates, the function MCMClogit is unable to
find the variable x1 (possibly x2,x3, and x4 also) in the data frame c.df.
Check the names of the variables in that data frame and make sure that the
names correspond to the formula specification.
Hope this helps,
Ravi.
----------------------------------------------------------------------------
-------
Ravi Varadhan, Ph.D.
Assistant Professor, The Center on Aging and Health
Division of Geriatric Medicine and Gerontology
Johns Hopkins University
Ph: (410) 502-2619
Fax: (410) 614-9625
Email: [EMAIL PROTECTED]
Webpage: http://www.jhsph.edu/agingandhealth/People/Faculty/Varadhan.html
----------------------------------------------------------------------------
--------
-----Original Message-----
From: [EMAIL PROTECTED]
[mailto:[EMAIL PROTECTED] On Behalf Of Anamika Chaudhuri
Sent: Friday, March 09, 2007 3:27 PM
To: [email protected]
Subject: [R] MCMC logit
Hi,
I have a dataset with the binary outcome Y(0,1) and 4 covariates
(X1,X@,X#,X$). I am trying to use MCMClogit to model logistic regression
using MCMC. I am getting an error where it doesnt identify the covariates
,although its reading in correctly. The dataset is a sample of actual
dataset. Below is my code:
> #######################
>
>
> #retreive data
> # considering four covariates
>
d.df=as.data.frame(read.table("c:/tina/phd/thesis/data/modified_data1.1.txt"
,header=T,sep=","))
> y=d.df[,ncol(d.df)]
> x=d.df[,1:4]
> c.df=cbind(y,x)
> #x=cbind(1,x)
> p <- ncol(c.df)
>
> # marginal log-prior of beta[]
> logpriorfun <- function(beta, mu, gshape, grate)
+ {
+ logprior = -p*log(2) + log(gamma(p+gshape)) - log(gamma(gshape))
+ + gshape*log(grate) - (p+gshape)* log(grate+sum(abs(beta)))
+ return(logprior)
+ }
> require(MCMCpack)
Loading required package: MCMCpack
Loading required package: coda
Loading required package: lattice
Loading required package: MASS
##
## Markov Chain Monte Carlo Package (MCMCpack)
## Copyright (C) 2003-2007 Andrew D. Martin and Kevin M. Quinn
##
## Support provided by the U.S. National Science Foundation
## (Grants SES-0350646 and SES-0350613)
##
[1] TRUE
Warning message:
package 'MASS' was built under R version 2.4.1
> a0 = 0.5
> b0 = 1
> mu0 = 0
> beta.init=list(c(0, rep(0.1,4)), c(0, rep(-0.1,4)), c(0, rep(0, 4)))
> burnin.cycles = 1000
> mcmc.cycles = 25000
> # three chains
> post.list <- lapply(beta.init, function(vec)
+ {
+ posterior <- MCMClogit(y~x1+x2+x3+x4, data=c.df, burnin=burnin.cycles,
mcmc=mcmc.cycles,
+ thin=5, tune=0.5, beta.start=vec, user.prior.density=logpriorfun,
logfun=TRUE,
+ mu=mu0, gshape=a0, grate=b0)
+ return(posterior)
+ })
Error in eval(expr, envir, enclos) : object "x1" not found
>
Any suggestions will be greatly appreciated.
Thanks,
Anamika
---------------------------------
We won't tell. Get more on shows you hate to love
[[alternative HTML version deleted]]
______________________________________________
[email protected] mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
---------------------------------
Bored stiff? Loosen up...
[[alternative HTML version deleted]]
______________________________________________
[email protected] mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.