EM.sclm: Censored Spatial Model Estimation via EM Algorithm

Description Usage Arguments Details Value Note Author(s) References See Also Examples

View source: R/EstEMspatial_USER.R

Description

This function returns the maximum likelihood (ML) estimates of the unknown parameters in Gaussian spatial models with censored/missing responses via the EM algorithm. It supports left, right, interval, or missing values in the dependent variable. It also computes the observed information matrix using the method developed by \insertCitelouis1982finding;textualRcppCensSpatial.

Usage

1
2
3
EM.sclm(y, x, cens, LI, LS, coords, init.phi, init.nugget,
  type = "exponential", kappa = 0, lower = c(0.01, 0.01), upper = c(30,
  30), MaxIter = 300, error = 1e-05, show.SE = TRUE)

Arguments

y

vector of responses.

x

design matrix.

cens

vector of censoring indicators. For each observation: 1 if censored/missing and 0 otherwise.

LI

lower limit of detection. For each observation: if non-censored =y, if left-censored/missing =-Inf, or =LOD if right/interval-censored.

LS

upper limit of detection. For each observation: if non-censored =y, if right-censored/missing =Inf, or =LOD if left/interval-censored.

coords

2D spatial coordinates.

init.phi

initial value for the spatial scaling parameter.

init.nugget

initial value for the nugget effect parameter.

type

type of spatial correlation function: 'exponential', 'gaussian', 'matern', and 'pow.exp' for exponential, gaussian, matern, and power exponential, respectively.

kappa

parameter for all spatial correlation functions. See CovMat.

lower, upper

vectors of lower and upper bounds for the optimization method. If unspecified, the default is c(0.01,0.01) for lower and c(30,30) for upper.

MaxIter

maximum number of iterations of the EM algorithm. By default =300.

error

maximum convergence error. By default =1e-5.

show.SE

TRUE or FALSE. It indicates if the standard errors should be estimated. By default =TRUE.

Details

The spatial Gaussian model is given by

Y = Xβ + ξ,

where Y is the n x 1 vector of response, X is the n x q design matrix, β is the q x 1 vector of regression coefficients to be estimated, and ξ is the error term which is normally distributed with zero-mean and covariance matrix Σ=σ^2 R(φ) + τ^2 I_n. We assume that Σ is non-singular and X has full rank \insertCitediggle2007springerRcppCensSpatial.

The estimation process was performed via the EM algorithm initially proposed by \insertCitedempster1977maximum;textualRcppCensSpatial. The conditional expectations are computed through function meanvarTMD available in package MomTrunc \insertCitegalarza2019momentsRcppCensSpatial.

Value

The function returns an object of class sclm which is a list given by:

Theta

estimated parameters in all iterations, θ = (β, σ^2, φ, τ^2).

theta

final estimation of θ = (β, σ^2, φ, τ^2).

beta

estimated β.

sigma2

estimated σ^2.

phi

estimated φ.

tau2

estimated τ^2.

EY

first moment for the truncated normal distribution computed in the last iteration.

EYY

second moment for the truncated normal distribution computed in the last iteration.

SE

vector of standard errors of θ = (β, σ^2, φ, τ^2).

InfMat

observed information matrix.

loglik

log-likelihood for the EM method.

AIC

Akaike information criterion.

BIC

Bayesian information criterion.

Iterations

number of iterations needed to converge.

ptime

processing time.

range

the effective range.

Note

The EM final estimates correspond to the estimates obtained at the last iteration of the EM algorithm.

To fit a regression model for non-censored data, just set cens as a vector of zeros.

Functions print, summary, and plot work for objects of class sclm.

Author(s)

Katherine L. Valeriano, Alejandro Ordonez, Christian E. Galarza and Larissa A. Matos.

References

\insertAllCited

See Also

SAEM.sclm, MCEM.sclm, predict.sclm

Examples

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
# Simulated example: 10% of left-censored observations
n = 50   # Test with another values for n
set.seed(1000)
coords = round(matrix(runif(2*n,0,15),n,2),5)
x = cbind(rnorm(n), runif(n))
data = rCensSp(c(-1,3),2,4,0.5,x,coords,"left",0.10,0,"gaussian",0)

fit = EM.sclm(y=data$yobs, x=data[,7:8], cens=data$cens, LI=data$LI,
             LS=data$LS, coords=data[,5:6], init.phi=3, init.nugget=1,
             type="gaussian", error=1e-4)
fit

RcppCensSpatial documentation built on Sept. 21, 2021, 5:07 p.m.