Error message running linear model of coregionalisation
145 views
Skip to first unread message
imi...@googlemail.com
Jun 16, 2014, 5:39:26 AM6/16/14
to r-inla-disc...@googlegroups.com
Hello,
I hope someone may be able to help me understand better the error message I am getting while trying to run a linear model of co-regionalisation. I'm not an expert and have only recently started learning INLA so it is possible that i have made a silly error!
I'm running Windows 7, R 3.0.1 (64 bit) with INLA version 0.0.1402298691 (I think I upgraded to the testing version, but the error was the same before doing this).
I'm trying to run a linear model of co-regionalisation (LMC). I have a dataset with 416 spatial locations (points). Each point has two binary outcome variables for several individuals (e.g. presence/absence of disease x and disease y, from population n). What I would like to do is to describe the shared spatial effect between disease x and disease y, using the LMC.
My code is below, followed by the output, including the error message.
I'd appreciate any input or advice regarding this. let me know if I should provide any additional information.
Thanks,
Nicola
#####
data<-read.csv("nicoladata.csv")
loc<-cbind(data$E, data$N)
plot(loc)
#create a mesh and A frame for the area
mesh <- inla.mesh.create.helper(points.domain=cbind(data$E, data$N),max.edge=c(0.04,50))
plot(mesh)
points(loc)
spde <- inla.spde2.matern(mesh, alpha=2)
Aj <- inla.spde.make.A(mesh=mesh, loc=loc)
#Outcome responses
X<-data$H_X
Y<-data$H_Y
#Create stack for disease X
stack.X<-inla.stack(data=list(X=X, y=cbind(X,NA)), A=list(Aj,1),effects=list(nodesX=1:spde$n.spde, intercept=rep(1,416)))
#Create stack for disease Y
stack.Y<-inla.stack(data=list(Y=Y, y=cbind(NA, Y)), A=list(Aj,1,1),effects=list(nodesY=1:spde$n.spde, randY=rep(1,416), intercept=rep(1,416)))
#Create stack for both disease X and disease Y
stack.both<-inla.stack(stack.X, stack.Y)
#Model disease X alone (this works ok)
formula.X<- X ~ intercept + f(nodesX, model = spde) - 1
result.X <- inla(formula.X,data=inla.stack.data(stack.X),family="binomial", Ntrials=data$H_X_N, control.predictor=list(A=inla.stack.A(stack.X)), control.compute=list(dic=TRUE), verbose=TRUE)
#Model disease Y alone (this works ok)
formula.Y<- Y ~ intercept + f(nodesY, model = spde) - 1
result.Y <- inla(formula.Y,data=inla.stack.data(stack.Y),family="binomial", Ntrials=data$H_Y_N, control.predictor=list(A=inla.stack.A(stack.Y)),control.compute=list(dic=TRUE), verbose=TRUE)
#Model both disease X and disease Y together (this does not work)
formulaboth <- y ~ intercept +
f(nodesX, model = spde) +
f(nodesY, copy="nodesX", hyper = list(beta=list(fixed=FALSE))) +
f(randY, model="iid", hyper=list(theta=list(prior="loggamma",param=c(1,0.01)))) - 1
resultboth <- inla(formulaboth,data=inla.stack.data(stack.both),family=c("binomial", "binomial"), control.compute=list(dic=TRUE), Ntrials=cbind(data$H_X_N, data$H_Y_N), control.predictor=list(A=inla.stack.A(stack.both)), verbose=TRUE)
############
############
OUTPUT FROM R:
hgid: a37e47c6b008 date: Mon Jun 09 08:53:00 2014 +0200
Report bugs to <he...@r-inla.org>
Processing file [C:/Users/NB1B06~1.SOT/AppData/Local/Temp/RtmpYDZL4r/file417c7d1c6f58/Model.ini] max_threads=[2]
inla_build...
number of sections=[9]
parse section=[0] name=[INLA.Model] type=[PROBLEM]
inla_parse_problem...
name=[INLA.Model]
openmp.strategy=[default]
store results in directory=[C:/Users/NB1B06~1.SOT/AppData/Local/Temp/RtmpYDZL4r/file417c7d1c6f58/results.files]
output:
cpo=[0]
po=[0]
dic=[1]
kld=[1]
mlik=[1]
q=[0]
graph=[0]
gdensity=[0]
hyperparameters=[1]
summary=[1]
return.marginals=[1]
nquantiles=[3] [ 0.025 0.5 0.975 ]
ncdf=[0] [ ]
parse section=[3] name=[Predictor] type=[PREDICTOR]
inla_parse_predictor ...
section=[Predictor]
dir=[predictor]
PRIOR->name=[loggamma]
PRIOR->from_theta=[function (x) <<NEWLINE>>exp(x)]
PRIOR->to_theta = [function (x) <<NEWLINE>>log(x)]
PRIOR->PARAMETERS=[1, 1e-005]
initialise log_precision[11]
fixed=[1]
user.scale=[1]
n=[774]
m=[832]
ndata=[832]
compute=[1]
read offsets from file=[C:/Users/NB1B06~1.SOT/AppData/Local/Temp/RtmpYDZL4r/file417c7d1c6f58/data.files/file417c7fdada1]
read n=[3212] entries from file=[C:/Users/NB1B06~1.SOT/AppData/Local/Temp/RtmpYDZL4r/file417c7d1c6f58/data.files/file417c7fdada1]
file=[C:/Users/NB1B06~1.SOT/AppData/Local/Temp/RtmpYDZL4r/file417c7d1c6f58/data.files/file417c7fdada1] 0/1606 (idx,y) = (0, 0)
file=[C:/Users/NB1B06~1.SOT/AppData/Local/Temp/RtmpYDZL4r/file417c7d1c6f58/data.files/file417c7fdada1] 1/1606 (idx,y) = (1, 0)
file=[C:/Users/NB1B06~1.SOT/AppData/Local/Temp/RtmpYDZL4r/file417c7d1c6f58/data.files/file417c7fdada1] 2/1606 (idx,y) = (2, 0)
file=[C:/Users/NB1B06~1.SOT/AppData/Local/Temp/RtmpYDZL4r/file417c7d1c6f58/data.files/file417c7fdada1] 3/1606 (idx,y) = (3, 0)
file=[C:/Users/NB1B06~1.SOT/AppData/Local/Temp/RtmpYDZL4r/file417c7d1c6f58/data.files/file417c7fdada1] 4/1606 (idx,y) = (4, 0)
file=[C:/Users/NB1B06~1.SOT/AppData/Local/Temp/RtmpYDZL4r/file417c7d1c6f58/data.files/file417c7fdada1] 5/1606 (idx,y) = (5, 0)
file=[C:/Users/NB1B06~1.SOT/AppData/Local/Temp/RtmpYDZL4r/file417c7d1c6f58/data.files/file417c7fdada1] 6/1606 (idx,y) = (6, 0)
file=[C:/Users/NB1B06~1.SOT/AppData/Local/Temp/RtmpYDZL4r/file417c7d1c6f58/data.files/file417c7fdada1] 7/1606 (idx,y) = (7, 0)
file=[C:/Users/NB1B06~1.SOT/AppData/Local/Temp/RtmpYDZL4r/file417c7d1c6f58/data.files/file417c7fdada1] 8/1606 (idx,y) = (8, 0)
file=[C:/Users/NB1B06~1.SOT/AppData/Local/Temp/RtmpYDZL4r/file417c7d1c6f58/data.files/file417c7fdada1] 9/1606 (idx,y) = (9, 0)
file=[C:/Users/NB1B06~1.SOT/AppData/Local/Temp/RtmpYDZL4r/file417c7d1c6f58/data.files/file417c7fdada1] 10/1606 (idx,y) = (10, 0)
file=[C:/Users/NB1B06~1.SOT/AppData/Local/Temp/RtmpYDZL4r/file417c7d1c6f58/data.files/file417c7fdada1] 11/1606 (idx,y) = (11, 0)
file=[C:/Users/NB1B06~1.SOT/AppData/Local/Temp/RtmpYDZL4r/file417c7d1c6f58/data.files/file417c7fdada1] 12/1606 (idx,y) = (12, 0)
file=[C:/Users/NB1B06~1.SOT/AppData/Local/Temp/RtmpYDZL4r/file417c7d1c6f58/data.files/file417c7fdada1] 13/1606 (idx,y) = (13, 0)
file=[C:/Users/NB1B06~1.SOT/AppData/Local/Temp/RtmpYDZL4r/file417c7d1c6f58/data.files/file417c7fdada1] 14/1606 (idx,y) = (14, 0)
file=[C:/Users/NB1B06~1.SOT/AppData/Local/Temp/RtmpYDZL4r/file417c7d1c6f58/data.files/file417c7fdada1] 15/1606 (idx,y) = (15, 0)
file=[C:/Users/NB1B06~1.SOT/AppData/Local/Temp/RtmpYDZL4r/file417c7d1c6f58/data.files/file417c7fdada1] 16/1606 (idx,y) = (16, 0)
file=[C:/Users/NB1B06~1.SOT/AppData/Local/Temp/RtmpYDZL4r/file417c7d1c6f58/data.files/file417c7fdada1] 17/1606 (idx,y) = (17, 0)
file=[C:/Users/NB1B06~1.SOT/AppData/Local/Temp/RtmpYDZL4r/file417c7d1c6f58/data.files/file417c7fdada1] 18/1606 (idx,y) = (18, 0)
file=[C:/Users/NB1B06~1.SOT/AppData/Local/Temp/RtmpYDZL4r/file417c7d1c6f58/data.files/file417c7fdada1] 19/1606 (idx,y) = (19, 0)
Aext=[C:/Users/NB1B06~1.SOT/AppData/Local/Temp/RtmpYDZL4r/file417c7d1c6f58/data.files/file417c701d24b0]
AextPrecision=[3.269e+006]
output:
summary=[1]
return.marginals=[1]
nquantiles=[3] [ 0.025 0.5 0.975 ]
ncdf=[0] [ ]
parse section=[1] name=[INLA.Data1] type=[DATA]
inla_parse_data [section 1]...
tag=[INLA.Data1]
family=[BINOMIAL]
likelihood=[BINOMIAL]
file->name=[C:/Users/NB1B06~1.SOT/AppData/Local/Temp/RtmpYDZL4r/file417c7d1c6f58/data.files/file417c6c564113]
file->name=[C:/Users/NB1B06~1.SOT/AppData/Local/Temp/RtmpYDZL4r/file417c7d1c6f58/data.files/file417c11b2848]
read n=[3328] entries from file=[C:/Users/NB1B06~1.SOT/AppData/Local/Temp/RtmpYDZL4r/file417c7d1c6f58/data.files/file417c6c564113]
*** ERROR *** inla_read_data_likelihood: file [C:/Users/NB1B06~1.SOT/AppData/Local/Temp/RtmpYDZL4r/file417c7d1c6f58/data.files/file417c6c564113] contains [3328] elements, which is not a multiple of [3]
Error in inla.inlaprogram.has.crashed() :
The inla-program exited with an error. Unless you interupted it yourself, please rerun with verbose=TRUE and check the output carefully.
If this does help; please contact the developers at <he...@r-inla.org>.
In addition: Warning message:
running command '"V:/R3.0.1_MNT/R/R-3.0.1/library/INLA/bin/windows/inla64.exe" -b -s -v "C:/Users/NB1B06~1.SOT/AppData/Local/Temp/RtmpYDZL4r/file417c7d1c6f58/Model.ini"' had status 1
>
I hope someone may be able to help me understand better the error message I am getting while trying to run a linear model of co-regionalisation. I'm not an expert and have only recently started learning INLA so it is possible that i have made a silly error!
I'm running Windows 7, R 3.0.1 (64 bit) with INLA version 0.0.1402298691 (I think I upgraded to the testing version, but the error was the same before doing this).
I'm trying to run a linear model of co-regionalisation (LMC). I have a dataset with 416 spatial locations (points). Each point has two binary outcome variables for several individuals (e.g. presence/absence of disease x and disease y, from population n). What I would like to do is to describe the shared spatial effect between disease x and disease y, using the LMC.
My code is below, followed by the output, including the error message.
I'd appreciate any input or advice regarding this. let me know if I should provide any additional information.
Thanks,
Nicola
#####
data<-read.csv("nicoladata.csv")
loc<-cbind(data$E, data$N)
plot(loc)
#create a mesh and A frame for the area
mesh <- inla.mesh.create.helper(points.domain=cbind(data$E, data$N),max.edge=c(0.04,50))
plot(mesh)
points(loc)
spde <- inla.spde2.matern(mesh, alpha=2)
Aj <- inla.spde.make.A(mesh=mesh, loc=loc)
#Outcome responses
X<-data$H_X
Y<-data$H_Y
#Create stack for disease X
stack.X<-inla.stack(data=list(X=X, y=cbind(X,NA)), A=list(Aj,1),effects=list(nodesX=1:spde$n.spde, intercept=rep(1,416)))
#Create stack for disease Y
stack.Y<-inla.stack(data=list(Y=Y, y=cbind(NA, Y)), A=list(Aj,1,1),effects=list(nodesY=1:spde$n.spde, randY=rep(1,416), intercept=rep(1,416)))
#Create stack for both disease X and disease Y
stack.both<-inla.stack(stack.X, stack.Y)
#Model disease X alone (this works ok)
formula.X<- X ~ intercept + f(nodesX, model = spde) - 1
result.X <- inla(formula.X,data=inla.stack.data(stack.X),family="binomial", Ntrials=data$H_X_N, control.predictor=list(A=inla.stack.A(stack.X)), control.compute=list(dic=TRUE), verbose=TRUE)
#Model disease Y alone (this works ok)
formula.Y<- Y ~ intercept + f(nodesY, model = spde) - 1
result.Y <- inla(formula.Y,data=inla.stack.data(stack.Y),family="binomial", Ntrials=data$H_Y_N, control.predictor=list(A=inla.stack.A(stack.Y)),control.compute=list(dic=TRUE), verbose=TRUE)
#Model both disease X and disease Y together (this does not work)
formulaboth <- y ~ intercept +
f(nodesX, model = spde) +
f(nodesY, copy="nodesX", hyper = list(beta=list(fixed=FALSE))) +
f(randY, model="iid", hyper=list(theta=list(prior="loggamma",param=c(1,0.01)))) - 1
resultboth <- inla(formulaboth,data=inla.stack.data(stack.both),family=c("binomial", "binomial"), control.compute=list(dic=TRUE), Ntrials=cbind(data$H_X_N, data$H_Y_N), control.predictor=list(A=inla.stack.A(stack.both)), verbose=TRUE)
############
############
OUTPUT FROM R:
hgid: a37e47c6b008 date: Mon Jun 09 08:53:00 2014 +0200
Report bugs to <he...@r-inla.org>
Processing file [C:/Users/NB1B06~1.SOT/AppData/Local/Temp/RtmpYDZL4r/file417c7d1c6f58/Model.ini] max_threads=[2]
inla_build...
number of sections=[9]
parse section=[0] name=[INLA.Model] type=[PROBLEM]
inla_parse_problem...
name=[INLA.Model]
openmp.strategy=[default]
store results in directory=[C:/Users/NB1B06~1.SOT/AppData/Local/Temp/RtmpYDZL4r/file417c7d1c6f58/results.files]
output:
cpo=[0]
po=[0]
dic=[1]
kld=[1]
mlik=[1]
q=[0]
graph=[0]
gdensity=[0]
hyperparameters=[1]
summary=[1]
return.marginals=[1]
nquantiles=[3] [ 0.025 0.5 0.975 ]
ncdf=[0] [ ]
parse section=[3] name=[Predictor] type=[PREDICTOR]
inla_parse_predictor ...
section=[Predictor]
dir=[predictor]
PRIOR->name=[loggamma]
PRIOR->from_theta=[function (x) <<NEWLINE>>exp(x)]
PRIOR->to_theta = [function (x) <<NEWLINE>>log(x)]
PRIOR->PARAMETERS=[1, 1e-005]
initialise log_precision[11]
fixed=[1]
user.scale=[1]
n=[774]
m=[832]
ndata=[832]
compute=[1]
read offsets from file=[C:/Users/NB1B06~1.SOT/AppData/Local/Temp/RtmpYDZL4r/file417c7d1c6f58/data.files/file417c7fdada1]
read n=[3212] entries from file=[C:/Users/NB1B06~1.SOT/AppData/Local/Temp/RtmpYDZL4r/file417c7d1c6f58/data.files/file417c7fdada1]
file=[C:/Users/NB1B06~1.SOT/AppData/Local/Temp/RtmpYDZL4r/file417c7d1c6f58/data.files/file417c7fdada1] 0/1606 (idx,y) = (0, 0)
file=[C:/Users/NB1B06~1.SOT/AppData/Local/Temp/RtmpYDZL4r/file417c7d1c6f58/data.files/file417c7fdada1] 1/1606 (idx,y) = (1, 0)
file=[C:/Users/NB1B06~1.SOT/AppData/Local/Temp/RtmpYDZL4r/file417c7d1c6f58/data.files/file417c7fdada1] 2/1606 (idx,y) = (2, 0)
file=[C:/Users/NB1B06~1.SOT/AppData/Local/Temp/RtmpYDZL4r/file417c7d1c6f58/data.files/file417c7fdada1] 3/1606 (idx,y) = (3, 0)
file=[C:/Users/NB1B06~1.SOT/AppData/Local/Temp/RtmpYDZL4r/file417c7d1c6f58/data.files/file417c7fdada1] 4/1606 (idx,y) = (4, 0)
file=[C:/Users/NB1B06~1.SOT/AppData/Local/Temp/RtmpYDZL4r/file417c7d1c6f58/data.files/file417c7fdada1] 5/1606 (idx,y) = (5, 0)
file=[C:/Users/NB1B06~1.SOT/AppData/Local/Temp/RtmpYDZL4r/file417c7d1c6f58/data.files/file417c7fdada1] 6/1606 (idx,y) = (6, 0)
file=[C:/Users/NB1B06~1.SOT/AppData/Local/Temp/RtmpYDZL4r/file417c7d1c6f58/data.files/file417c7fdada1] 7/1606 (idx,y) = (7, 0)
file=[C:/Users/NB1B06~1.SOT/AppData/Local/Temp/RtmpYDZL4r/file417c7d1c6f58/data.files/file417c7fdada1] 8/1606 (idx,y) = (8, 0)
file=[C:/Users/NB1B06~1.SOT/AppData/Local/Temp/RtmpYDZL4r/file417c7d1c6f58/data.files/file417c7fdada1] 9/1606 (idx,y) = (9, 0)
file=[C:/Users/NB1B06~1.SOT/AppData/Local/Temp/RtmpYDZL4r/file417c7d1c6f58/data.files/file417c7fdada1] 10/1606 (idx,y) = (10, 0)
file=[C:/Users/NB1B06~1.SOT/AppData/Local/Temp/RtmpYDZL4r/file417c7d1c6f58/data.files/file417c7fdada1] 11/1606 (idx,y) = (11, 0)
file=[C:/Users/NB1B06~1.SOT/AppData/Local/Temp/RtmpYDZL4r/file417c7d1c6f58/data.files/file417c7fdada1] 12/1606 (idx,y) = (12, 0)
file=[C:/Users/NB1B06~1.SOT/AppData/Local/Temp/RtmpYDZL4r/file417c7d1c6f58/data.files/file417c7fdada1] 13/1606 (idx,y) = (13, 0)
file=[C:/Users/NB1B06~1.SOT/AppData/Local/Temp/RtmpYDZL4r/file417c7d1c6f58/data.files/file417c7fdada1] 14/1606 (idx,y) = (14, 0)
file=[C:/Users/NB1B06~1.SOT/AppData/Local/Temp/RtmpYDZL4r/file417c7d1c6f58/data.files/file417c7fdada1] 15/1606 (idx,y) = (15, 0)
file=[C:/Users/NB1B06~1.SOT/AppData/Local/Temp/RtmpYDZL4r/file417c7d1c6f58/data.files/file417c7fdada1] 16/1606 (idx,y) = (16, 0)
file=[C:/Users/NB1B06~1.SOT/AppData/Local/Temp/RtmpYDZL4r/file417c7d1c6f58/data.files/file417c7fdada1] 17/1606 (idx,y) = (17, 0)
file=[C:/Users/NB1B06~1.SOT/AppData/Local/Temp/RtmpYDZL4r/file417c7d1c6f58/data.files/file417c7fdada1] 18/1606 (idx,y) = (18, 0)
file=[C:/Users/NB1B06~1.SOT/AppData/Local/Temp/RtmpYDZL4r/file417c7d1c6f58/data.files/file417c7fdada1] 19/1606 (idx,y) = (19, 0)
Aext=[C:/Users/NB1B06~1.SOT/AppData/Local/Temp/RtmpYDZL4r/file417c7d1c6f58/data.files/file417c701d24b0]
AextPrecision=[3.269e+006]
output:
summary=[1]
return.marginals=[1]
nquantiles=[3] [ 0.025 0.5 0.975 ]
ncdf=[0] [ ]
parse section=[1] name=[INLA.Data1] type=[DATA]
inla_parse_data [section 1]...
tag=[INLA.Data1]
family=[BINOMIAL]
likelihood=[BINOMIAL]
file->name=[C:/Users/NB1B06~1.SOT/AppData/Local/Temp/RtmpYDZL4r/file417c7d1c6f58/data.files/file417c6c564113]
file->name=[C:/Users/NB1B06~1.SOT/AppData/Local/Temp/RtmpYDZL4r/file417c7d1c6f58/data.files/file417c11b2848]
read n=[3328] entries from file=[C:/Users/NB1B06~1.SOT/AppData/Local/Temp/RtmpYDZL4r/file417c7d1c6f58/data.files/file417c6c564113]
*** ERROR *** inla_read_data_likelihood: file [C:/Users/NB1B06~1.SOT/AppData/Local/Temp/RtmpYDZL4r/file417c7d1c6f58/data.files/file417c6c564113] contains [3328] elements, which is not a multiple of [3]
Error in inla.inlaprogram.has.crashed() :
The inla-program exited with an error. Unless you interupted it yourself, please rerun with verbose=TRUE and check the output carefully.
If this does help; please contact the developers at <he...@r-inla.org>.
In addition: Warning message:
running command '"V:/R3.0.1_MNT/R/R-3.0.1/library/INLA/bin/windows/inla64.exe" -b -s -v "C:/Users/NB1B06~1.SOT/AppData/Local/Temp/RtmpYDZL4r/file417c7d1c6f58/Model.ini"' had status 1
>
INLA help
Jun 16, 2014, 5:44:15 AM6/16/14
to imi...@googlemail.com, r-inla-disc...@googlegroups.com
On Mon, 2014-06-16 at 02:39 -0700, imi...@googlemail.com wrote:
> data<-read.csv("nicoladata.csv")
Hi,
can yu send this csv-file as well (or a coded version of it with same
dimension) ?
H
--
Håvard Rue
he...@r-inla.org
> data<-read.csv("nicoladata.csv")
Hi,
can yu send this csv-file as well (or a coded version of it with same
dimension) ?
H
--
Håvard Rue
he...@r-inla.org
Finn Lindgren
Jun 16, 2014, 5:57:34 AM6/16/14
to r-inla-disc...@googlegroups.com
On 16/06/14 10:39, imi...@googlemail.com wrote:
> resultboth <-
> inla(formulaboth,data=inla.stack.data(stack.both),family=c("binomial",
> "binomial"), control.compute=list(dic=TRUE), Ntrials=cbind(data$H_X_N,
> data$H_Y_N), control.predictor=list(A=inla.stack.A(stack.both)),
> verbose=TRUE)
I believe your Ntrials specification here is wrong, as it has fewer rows
> resultboth <-
> inla(formulaboth,data=inla.stack.data(stack.both),family=c("binomial",
> "binomial"), control.compute=list(dic=TRUE), Ntrials=cbind(data$H_X_N,
> data$H_Y_N), control.predictor=list(A=inla.stack.A(stack.both)),
> verbose=TRUE)
that the combined data. In this case, I think the correct solutions is
to use the vector (not matrix)
Ntrials = c(data$H_X_N, data$H_Y_N)
In general, I recommend keeping track of Ntrials in the "data" section
of the stacks to make sure it gets the correct size, just as you do for
the data itself.
Finn
Message has been deleted
imi...@googlemail.com
Jun 16, 2014, 6:33:34 AM6/16/14
to r-inla-disc...@googlegroups.com
Fantastic - that has indeed solved the problem.
Thanks very much for your input!
Nicola
Reply all
Reply to author
Forward
0 new messages