I'm having trouble with adapt. I'm trying to use it in a Bayesian setting,
to integrate the posterior distribution, and to find posterior means. I
tried using the following script, and things went ok:
data = rnorm(100,0.2,1.1)
data = c(data,rnorm(10,3,1))
data = data[abs(data)<2*sd(data)]
prior = function(x){
dgamma(x[2],shape=2,scale=1)*dnorm(x[1],0,.5)
}
liklihood = function(x,val){
prod(dnorm(val,m=x[1],sd=sqrt(.5/x[2])))#/(pnorm(2,x[1],x[2])-pnorm(-2,x[1],x[2]))^length(val)
}
unsc.post = function(x,val){
prior(x)*liklihood(x,val)
}
cons = adapt(2,c(-5,0),c(5,10),f=unsc.post,min=1e+04,max=5e+05,val=data)
However, if I try to add a skewness parameter using sn as follows, I get an
error.
data = rnorm(20,0.2,1.1)
data = c(data,rnorm(1,3,1))
data = data[abs(data)<2*sd(data)]
prior = function(x){
dgamma(x[2],shape=2,scale=1)*dnorm(x[1],0,.5)*dnorm(x[3],0,.1)
}
liklihood = function(x,val){
prod(dsn(val,loc=x[1],sc=sqrt(.5/x[2]),sh=x[3]))
}
unsc.post = function(x,val){
prior(x)*liklihood(x,val)
}
cons = adapt(2,c(-5,0,-.5),c(5,10,.5),f=unsc.post,max=5e+05,val=data)
This produces the error:
*** caught segfault ***
address 0x7eaa8580, cause 'memory not mapped'
Traceback:
1: .C("cadapt", as.integer(ndim), as.double(lower), as.double(upper),
minpts = as.integer(minpts), maxpts = as.integer(maxpts), ff, rho =
environment(), as.double(eps), relerr = double(1), lenwrk =
as.integer(lenwrk),
value = double(1), ifail = integer(1), PACKAGE = "adapt")
2: adapt(2, c(-5, 0, -0.5), c(5, 10, 0.5), f = unsc.post, max = 5e+05,
val = data)
On a perhaps related note, if I increase the length of data to 1000+10 in
the original example, I get a garbage answer and:
Warning message:
Ifail=2, lenwrk was too small. -- fix adapt() !
Check the returned relerr! in: adapt(2, c(-5, 0), c(5, 10), f = unsc.post,
min = 10000, max = 5e+05.
And If I uncomment the line there with data size 100+10, adapt takes forever
and gives me an error saying maxpts is too low (I don't think this is
related to the previous problem, but I thought I'd provide the information
in case it is useful).
Any help would be greatly appreciated.
Thanks,
Omkar Muralidharan
[[alternative HTML version deleted]]
______________________________________________
R-help@r-project.org 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.
