FirstOrderCCL: Semi-parametric multinomial logistic regression for...

In kristianhessellund/Multilogreg: Semi-parametric models for multivariate point pattern data

Description

Semi-parametric multinomial logistic regression for multivariate point pattern data.

Usage

FirstOrderCCL(X, Beta0, covariate)

Arguments

X	A multivariate point pattern.
Beta0	Initial value for Beta.
covariate	A matrix or image of covariates.The rows corresponds to the points in the point pattern and the columns indicates the corresponding covariates that are observed at the location of the point pattern.

Value

Return estimated first order parameters of X

Author(s)

Kristian Bjørn Hessellund, Ganggang Xu, Yongtao Guan and Rasmus Waagepetersen.

References

Hessellund, K. B., Xu, G., Guan, Y. and Waagepetersen, R. (2020) Semi-parametric multinomial logistic regression for multivariate point pattern data.

Examples

### First order analysis of Washington D.C. street crimes ###

nspecies <- length(unique(markedprocess$marks))
ncovs <- dim(Z_process)[2]

# Initial parameters for street crime data
Beta0=matrix(runif(nspecies*ncovs,-0.5,0.5),nrow =nspecies)

# Choosing the last type of street crime as the baseline
Beta0 <- as.matrix(t(scale(Beta0,center=Beta0[nspecies,],scale=FALSE)))

# Parameter estimation
Betahat=FirstOrderCCL(X=markedprocess,Beta0 = Beta0,covariate = Z_process)

# Parameter estimates
colnames(Betahat$betahat)=c("Burglary", "Assault w. weapon", "Motor v. theft",
"Theft f. auto", "Robbery", "Other theft (baseline)")
rownames(Betahat$betahat)=c("Intercept", "%African", "%Hispanic", "%Males",
 "Median income", "%Household","%Bachelor", "Distance to police station")
Betahat

### Non-parametric estimates of cross PCF ratios ###

rseq = seq(0,3000,length=100)
cpcf=CrossPCF(X = markedprocess,covariate = Z_process,Beta = Betahat$betahat,r = rseq
,bwd = 200,cut = NULL)

plot(rseq,rep(1,100),type="l",col=1,ylim=c(0,2),xlab="r",ylab="",main="PCFs")
for(i in 1:(nspecies-1)){
  lines(rseq,cpcf$cPCF.m[i,i,],,col=i)
}

plot(rseq,rep(1,100),type="l",col=1,ylim=c(0,2),xlab="r",ylab="",main="PCFs")
for(i in 1:(nspecies-1)){
  for(j in 1:(nspecies-1)){
    if(i==j){next}
    lines(rseq,cpcf$cPCF.m[i,j,],,col=1)
  }
}

### Standard errors for Beta ###

var.est=Firstorder_var(X=markedprocess,covariate = Z_process,Beta=Betahat$betahat,
r = seq(0,3000,length=100),bwd = 200,cPCF = cpcf$cPCF.m,method = "pool")

std.error=matrix(sqrt(diag(var.est$COV)),nrow=ncovs)
colnames(std.error)=c("Burglary", "Assault w. weapon", "Motor v. theft","Theft f. auto", "Robbery")
rownames(std.error)=c("Intercept", "%African", "%Hispanic", "%Males",
 "Median income", "%Household","%Bachelor", "Distance to police station")
std.error

###################################################################

### First order analysis of simulated LGCP ###

# Size of the observation window
n.x <- n.y <- 1
xx=seq(0,n.x,length=100)
yy=seq(0,n.x,length=100)

# Simulating a covariate
cov <- as.matrix(RFsimulate(RMexp(var=1,scale=0.05), x=xx, y=yy, grid=TRUE))

# Simulating the background intensity
gamma <- 0.5
background <- as.matrix(RFsimulate(RMgauss(var=1,scale=0.2), x=xx, y=yy, grid=TRUE))*gamma-gamma^2/2

#Set up parameters
beta1 <- c(0.1,0.2,0.3,0.4,0.5)
beta2 <- c(-0.1,-0.2,0,0.1,0.2)
beta2 <- as.matrix(beta2)

# Parameters in the LGCP

alpha <- matrix(c(0.5,-1,0.5,0,-1,0,0,0.5,0,0.5),nrow=5,byrow=TRUE)
xi    <- c(0.02,0.03)
sigma <- matrix(c(sqrt(0.5),sqrt(0.5),sqrt(0.5),sqrt(0.5),sqrt(0.5)),ncol=1)
phi   <- matrix(c(0.02,0.02,0.03,0.03,0.04),ncol=1)

n.window <- n.x
n.points <- c(400,400,400,400,400)

# Simulation of a multivariate LGCP
X <- sim_lgcp_multi(basecov=background,covariate=cov,betas=beta2,alphas=alpha,xi=xi,
sigma=sigma,phis=phi, n.window=n.window,n.points=n.points,beta0s=beta1)

nspecies  <- dim(beta2)[1]
Beta.true <- cbind(beta1,beta2)
Beta.true <- as.matrix(t(scale(Beta.true,center=Beta.true[nspecies,],scale=FALSE)))
Betahat   <- FirstOrderCCL(X=X$markedprocess,Beta0=Beta.true,covariate = as.im(t(cov)))
Betahat

kristianhessellund/Multilogreg documentation built on Jan. 1, 2021, 7:23 a.m.