Problem with singular matrix errosarlm

Hi,

I am trying to fit a spatial regression model with errorsarlm in spdep, 
but I have been running into some problems.

I can get weights style W and B to work with method LU however not S or 
C I get the following error (traceback below): Error in 
solve.default(-(mat), tol.solve = tol.solve) : system is computationally 
singular: reciprocal condition number = 0
My guess would be that the interval= argument should be given for S and C. 
In the LU case, the interval is set to -1,0.999 internally, and for each 
of these cases, the spatial coefficient may be being returned outside its 
domain, thus leading to failure in inverting the matrix returned by the 
numerical Hessian. Something similar may happen with Matrix - the Cholesky 
method.
I have rescaled the RHS variables to mean 0 and standard deviation 1. I 
have also tried removing all but two of the RHS variables which are 
definitely not collinear to see if that was the problem but I get the 
same error.

Each data point represents a 1by1 degree cell for the world's oceans, 
because I am using dnearneigh points at the equator have far fewer 
neighbours than those nearer the poles. Could this unequal number of 
neighbours be the reason for the error?

Also this model has 4647 points and takes around 50minutes to run which 
seems like a long time do you think this is due to the number of 
neighbours? The full data set has 31609 so I imagine this is going to 
take an extremely long time.

This is almost certainly caused by your choice of distance weights - you 
have 52% dense matrices, so lots of operations are being done on almost 
dense matrices. With 4647 points, you can use the eigen method too, which 
would give the domain interval exactly. You seem to have a mean number of 
links of 2461, so something rather odd is going on. The 1by1 degree cell 
for the world's land areas used in Bivand, Hauke & Kossowski (Geographical 
Analysis, 2013) with about 15K cells doesn't cause problems, but is very 
sparse.

Hope this clarifies,

Roger
Code:

chl<-(chl-mean(chl,na.rm=T))/sd(chl)
sst<-(sst-mean(sst,na.rm=T))/sd(sst)
chlm<-(chlm-mean(chlm,na.rm=T))/sd(chlm)
sstm<-(sstm-mean(sstm,na.rm=T))/sd(sstm)
depth<-(depth-mean(depth,na.rm=T))/sd(depth)
lats<-(lats-mean(lats,na.rm=T))/sd(lats)
longs<-(longs-mean(longs,na.rm=T))/sd(longs)

ols<-lm(chl~sst+depth+chlm+sstm+lats+longs)

coords<-cbind(longs,lats)
coords<-as.matrix(coords)

#Define neighbourhood here max is 200km
nn<-dnearneigh(coords,0,200,longlat = TRUE)

nn
Neighbour list object:
Number of regions: 4647
Number of nonzero links: 11438816
Percentage nonzero weights: 52.9707
Average number of links: 2461.549

#get inverse distances
dsts <- nbdists(nn, coords)
idw <- lapply(dsts, function(x) 1/(x))

#Spatial weights
nnweights<-nb2listw(nn, glist=idw, style='S')

sem0_500<-errorsarlm(formula(ols), listw=nnweights,na.action=na.omit, method='LU',
+               tol.solve=1e-16, control=list(returnHcov=FALSE))
Error in solve.default(-(mat), tol.solve = tol.solve) :
 system is computationally singular: reciprocal condition number = 0

traceback()
5: solve.default(-(mat), tol.solve = tol.solve)
4: solve(-(mat), tol.solve = tol.solve)
3: solve(-(mat), tol.solve = tol.solve)
2: getVmate(coefs, env, s2, trs, tol.solve = tol.solve, optim = con$optimHess,
      optimM = con$optimHessMethod)
1: errorsarlm(formula(ols), listw = nnweights, na.action = na.omit,
      method = "LU", tol.solve = 1e-16, control = list(returnHcov = FALSE))

when I try to run it with different method I got a different error

sem0_500<-errorsarlm(formula(ols), listw=nnweights,na.action=na.omit, method='Matrix',
+  tol.solve=1e-16, control=list(returnHcov=FALSE))
Warning message: In .local(object, ...) :
 Cholmod warning 'matrix not positive definite' at file ../Supernodal/t_cholmod_super_numeric.c, line 729

Error in determinant(update(nChol, nW, 1/coef)) :
 error in evaluating the argument 'x' in selecting a method for function 'determinant': Error in .local(object, ...) : CHOLMOD factorization was unsuccessful

traceback()
8: determinant(update(nChol, nW, 1/coef))
7: ifelse(coef > b, n * log(coef) + (.f * c(determinant(update(nChol,
      nW, 1/coef))$modulus)), ifelse(coef < a, n * log(-(coef)) +
      (.f * c(determinant(update(pChol, csrw, 1/(-coef)))$modulus)),
      0))
6: Matrix_ldet(coef, env, which = which)
5: do_ldet(lambda, env)
4: f(arg, ...)
3: (function (arg)
  -f(arg, ...))(0.930150734356483)
2: optimize(sar.error.f, interval = interval, maximum = TRUE, tol = con$tol.opt,
      env = env)
1: errorsarlm(formula(ols), listw = nnweights, na.action = na.omit,
      method = "Matrix", tol.solve = 1e-16, control = list(returnHcov = FALSE))

sessionInfo()
R version 3.0.2 (2013-09-25)
Platform: x86_64-pc-linux-gnu (64-bit)

locale:
[1] LC_CTYPE=en_GB.UTF-8       LC_NUMERIC=C
[3] LC_TIME=en_GB.UTF-8        LC_COLLATE=en_GB.UTF-8
[5] LC_MONETARY=en_GB.UTF-8    LC_MESSAGES=en_GB.UTF-8
[7] LC_PAPER=en_GB.UTF-8       LC_NAME=C
[9] LC_ADDRESS=C               LC_TELEPHONE=C
[11] LC_MEASUREMENT=en_GB.UTF-8 LC_IDENTIFICATION=C

attached base packages:
[1] stats     graphics  grDevices utils     datasets  methods   base

other attached packages:
[1] ncf_1.1-4    spdep_0.5-68 Matrix_1.1-0 sp_1.0-14

loaded via a namespace (and not attached):
[1] boot_1.3-9      coda_0.16-1     deldir_0.1-1    grid_3.0.2
[5] lattice_0.20-24 LearnBayes_2.12 MASS_7.3-29     nlme_3.1-111
[9] splines_3.0.2

I feel I am doing something wrong when specifying the weights matrix but I can't work out what. Any help would be hugely appreciated.

Thanks,

Emma

PhD student
Imperial College London

_______________________________________________
R-sig-Geo mailing list
R-sig-Geo at r-project.org
https://stat.ethz.ch/mailman/listinfo/r-sig-geo

Roger Bivand
Department of Economics, Norwegian School of Economics,
Helleveien 30, N-5045 Bergen, Norway.
voice: +47 55 95 93 55; fax +47 55 95 95 43
e-mail: Roger.Bivand at nhh.no

Problem with singular matrix errosarlm

Thread (3 messages)