fsglmm.discrete: Fitting Projection Based Laplace Approximation for Spatial...
In fastLaplace: A Fast Laplace Method for Spatial Generalized Linear Mixed Model

Description Usage Arguments Value Examples

View source: R/fll_discrete.R

fsglmm.discrete

fsglmm.discrete(
  formula,
  inits,
  data,
  family,
  ntrial = 1,
  method.optim,
  method.integrate,
  rank = NULL,
  A,
  offset = NA
)

`formula`	an object of class "formula."
`inits`	starting values for the parameters.
`data`	a data frame containing variables in the model.
`family`	a character string of the error distribution and link function to be used in the model.
`ntrial`	a numeric vector for binomial model.
`method.optim`	the method to be used for outer optimization. "CG" for Conjugate Gradient Method.
`method.integrate`	the method to be used for inner optimization. "NR" for Newton Raphson Method.
`rank`	an integer of 'rank' to be used for projections. Default is 5 percent of observations.
`A`	an adjacency matrix
`offset`	this is used to specify an a priori a known component to be included in the linear predictor during fitting.

a list containing the following components:

summary a summary of the fitted model

mle2 an object of class "mle2"

Delta a matrix containing the estimated random effects of the reduced dimensional model.

M the projection matrix used.

if(requireNamespace("ngspatial")&
requireNamespace("mgcv")){
n = 30
A = ngspatial::adjacency.matrix(n)
Q = diag(rowSums(A),n^2) - A
x = rep(0:(n - 1) / (n - 1), times = n)
y = rep(0:(n - 1) / (n - 1), each = n)
X = cbind(x, y)
beta = c(1, 1)
P.perp = diag(1,n^2) - X%*%solve(t(X)%*%X)%*%t(X)
eig = eigen(P.perp %*% A %*% P.perp)
eigenvalues = eig$values
q = 400
M = eig$vectors[,c(1:q)]
Q.s = t(M) %*% Q %*% M
tau = 6
Sigma = solve(tau*Q.s)
set.seed(1)
delta.s = mgcv::rmvn(1, rep(0,q), Sigma)
lambda = exp( X%*%beta + M%*%delta.s )
Z = c()
for(j in 1:n^2){Z[j] = rpois(1,lambda[j])}
Y = as.matrix(Z,ncol=1)
data = data.frame("Y"=Y,"X"=X)
colnames(data) = c("Y","X1","X2")
linmod <- glm(Y~-1+X1+X2,data=data,family="poisson") # Find starting values
linmod$coefficients
starting <- c(linmod$coefficients,"logtau"=log(1/var(linmod$residuals)) )

result.pois.disc <- fsglmm.discrete(Y~-1+X1+X2, inits = starting, data=data,
family="poisson",ntrial=1, method.optim="BFGS", method.integrate="NR",rank=50, A=A)
}