CAR model with generic1 or rgeneric?

[Stefano Castruccio]

Good evening RINLA community, 

                 I have been recently trying to use a CAR model for areal data according to the parametrization of the precision matrix as 1/tau*(I-rho*W), with rho in [0,1). Now, since the default parametrization for CAR is "besagproper", which is slightly different, one option is to build a precision matrix from scratch with rgeneric (which as been in fact been done here: https://becarioprecario.bitbucket.io/inla-gitbook/ch-newmodels.html#), but the generic1 option seems to me a much simpler option, since it deals with precision matrices of the form 1/tau*(I-beta/lambda_max*C), with beta in [0,1) and lambda_max maximum eigenvalue of C. So wouldn't implementing it be as simple as specifying generic1 with Cmatrix=W*lambda_max? Since I could not find this explicitly written anywhere, so it's either too obvious or I am missing something. 

Another related question is on SAR models: is there a simpler way than hard coding the precision matrix with rgeneric?

[Helpdesk]

you can do that, yes.

I see lambda.max is computed internally in the C-program, kind of sub-
optimal, through computing all eigenvalues. its possible to use like
'RSpectre' to compute directly only the max eigenvalue, which might be
better for larger models.

If that becomes an issue, then let me know and I can redo the code.

H

> --
> You received this message because you are subscribed to the Google
> Groups "R-inla discussion group" group.
> To unsubscribe from this group and stop receiving emails from it, send
> an email to r-inla-discussion...@googlegroups.com.
> To view this discussion on the web, visit
> https://groups.google.com/d/msgid/r-inla-discussion-group/68d54c14-77c3-417a-8b06-22084e78e191n%40googlegroups.com
> .

--
Håvard Rue
he...@r-inla.org

[Finn Lindgren]

I wrote a long answer and then accidentally deleted it...

Below is a recreation of it.

SAR(2) models have exactly the same precision structure as Matern SPDE models with alpha=2 (2nd order PDE operator, 4th order precision operator) and can be implemented using the generic structure defined in the attached "spde2" model definition implemented in inla.spde2.generic().

Example: For a SAR(2) with precision structure Q = tau * (I-beta*W)'*(I-beta*W) (with tau as a precision scaling parameter), the only thing that's limited in inla.spde2.generic is the allowed domain for the beta parameter.

Rewriting the precision as

  Q = tau * beta^2 * (beta^-2 * I - beta * (W' +W) + W'*W)

and comparing to the spde2 model structure with transform="log" for the phi2 parameter,

  Q = e^(2*phi0) * (e^(2*phi1) * M0 + (2*exp(phi2)-1) * exp(phi1) * (M1'+M1) + M2)

gives either (if my math is right)

  phi0 = log(tau*beta^2) / 2 = log(tau) / 2 + log(beta) = theta_tau/2 + theta_beta

  phi1 = -log(beta) = -theta_beta

  phi2 = -Inf

  M0 = I

  M1 = W

  M2 = W'*W

  B0 = cbind(0, 1/2, 1)

  B1 = cbind(0, 0, -1) and

  B2 = cbind(-Inf, 0, 0)

or

  phi0 = log(tau*beta^2) / 2 = log(tau) / 2 + log(beta) = theta_tau/2 + theta_beta

  phi1 = -log(beta) = -theta_beta

  phi2 = 0

  M0 = I

  M1 = -W

  M2 = W'*W

  B0 = cbind(0, 1/2, 1)

  B1 = cbind(0, 0, -1)

  B2 = cbind(0, 0, 0)

Both versions give theta_tau = log(tau) and theta_beta = log(beta) as internal INLA paramters.

  

The -Inf in the first version would need to be replaced by a finite value, so the second option looks more reasonable.

But I haven't considered exactly what the parameterisation domain can be with other choices; the versions above both imply beta > 0,

but for discrete SAR(1), all beta values are valid (but you get intrinsic fields when 1/beta=eigenvalue of W).

For beta < 0, just change the sign of W...

Finn

--
You received this message because you are subscribed to the Google Groups "R-inla discussion group" group.
To unsubscribe from this group and stop receiving emails from it, send an email to r-inla-discussion...@googlegroups.com.
To view this discussion on the web, visit https://groups.google.com/d/msgid/r-inla-discussion-group/68d54c14-77c3-417a-8b06-22084e78e191n%40googlegroups.com.

--

Finn Lindgren
email: finn.l...@gmail.com

[Stefano Castruccio]

Thanks Håvard, 

            for my case it's not a big issue as I only have about a hundred areas.