[R-sig-Geo] GSTAT Co-Kriging: How can i fix a non positive definte correlation matrix

Wed Jan 7 09:08:57 CET 2015

First: you don't provide us with a reproducible example. I know this is
a bit of work for your case, but it may be worth it. For instance, in
the text below several R commands appear on a single line, which is
syntactically an error.

Second: the error message you get may come from different reasons, one
of them being near-perfect correlation among variables, but you seem to
have prevented this with argument correct.diagonal; have you tried a
slightly larger value? Another cause may be duplicate or near duplicate
observations; sp::zerodist may be useful here.

On 01/06/2015 10:52 PM, Roberts, Blake wrote:
> QUESTION: Can anyone give me advice on how to resolve an error message regarding a non positive definite correlation matrix so that the computed negative cross-variogram can be used to perform co-kriging?
> 
> I am performing co-kriging on negatively correlated parameters. 
> 
> data = dtrd 
> 
> x	y	n.ts	n.ae
> 3.2	6.4	-0.052514693	0.055104755
> 6.4	6.4	-0.375234492	-0.861381682
> 9.6	6.4	1.362810557	-0.23577909
> 12.8	6.4	-0.333087778	0.134741672
> 16	6.4	0.699916885	-2.348188189
> 22.4	6.4	0.605552409	-0.75492334
> 25.6	6.4	0.406868496	0.329847607
> 3.2	9.6	-0.529754386	0.251549504
> 6.4	9.6	-0.416509276	1.094807283
> 9.6	9.6	-0.450700276	0.988002955
> 12.8	9.6	-0.732937251	1.344413512
> 16	9.6	0.377457633	-0.542187169
> 19.2	9.6	-0.989029229	2.002433474
> 22.4	9.6	-1.259551904	-0.895859924
> 25.6	9.6	-0.470220669	0.370809323
> 3.2	12.8	1.046426438	0.168121785
> 6.4	12.8	-1.440639861	1.423433672
> 9.6	12.8	-0.899563555	0.83844681
> 12.8	12.8	-1.287567164	-0.086478673
> 16	12.8	-0.816931437	1.156633451
> 19.2	12.8	-0.782655308	0.299641538
> 22.4	12.8	0.048342822	-0.763091614
> 3.2	16	-0.5581308	-0.130509134
> 6.4	16	-0.217874377	0.815852429
> 9.6	16	-0.272702769	0.906621066
> 12.8	16	-0.996064	0.274459785
> 16	16	-1.97520507	2.017626766
> 19.2	16	-0.718163458	0.444173367
> 22.4	16	0.337266312	0.831151021
> 25.6	16	-0.565668608	0.513114809
> 3.2	19.2	0.651736313	-0.055114304
> 6.4	19.2	-0.715751462	0.205488257
> 9.6	19.2	0.650820529	0.558032798
> 12.8	19.2	-1.184380204	0.612006591
> 16	19.2	0.132388601	0.462059331
> 19.2	19.2	1.08502868	0.399240335
> 22.4	19.2	-0.85196341	0.113079949
> 25.6	19.2	0.099887288	-1.922971577
> 3.2	22.4	-0.233091874	-0.287556625
> 6.4	22.4	-0.989196908	0.10047403
> 9.6	22.4	-0.025054063	0.384754057
> 12.8	22.4	-0.394999955	-1.272051846
> 16	22.4	0.049839033	-2.568551091
> 19.2	22.4	1.327518037	-0.974203156
> 22.4	22.4	-0.421157849	-0.0101211
> 25.6	22.4	-0.178717505	0.408401438
> 3.2	25.6	-0.065413063	-0.195765074
> 6.4	25.6	1.167410569	-1.626236059
> 9.6	25.6	0.428968564	-1.287922066
> 12.8	25.6	1.550074254	-1.011990828
> 16	25.6	2.931334909	-1.872352342
> 19.2	25.6	4.126464345	-2.254471152
> 22.4	25.6	0.44561262	-0.709554064
> 25.6	25.6	0.543591219	1.03999109
> 3.2	28.8	0.716228164	-0.671315106
> 6.4	28.8	1.856183532	-1.403466298
> 9.6	28.8	0.774644882	-0.952601604
> 12.8	28.8	0.425336383	-0.776073605
> 16	28.8	0.294908064	-0.00019281
> 19.2	28.8	2.38260471	-0.778811406
> 22.4	28.8	-0.315391214	-0.223383771
> 25.6	28.8	1.136916243	-0.238938091
> 3.2	32	-0.233091874	-0.520269718
> 6.4	32	-0.352017426	-0.556041646
> 9.6	32	-0.933099317	-0.441354858
> 16	32	-0.865945243	-1.011088256
> 19.2	32	-0.594339105	-0.399550741
> 22.4	32	-0.653328511	1.213119402
> 25.6	32	0.076670222	0.819568017
> 3.2	35.2	-0.14022361	0.166527242
> 6.4	35.2	-0.545492977	-0.613535469
> 9.6	35.2	-0.084386566	-0.423754708
> 12.8	35.2	2.107283841	-0.891753223
> 16	35.2	-1.645006796	0.566968259
> 22.4	35.2	0.512684145	0.515250895
> 25.6	35.2	0.517794479	-1.144578624
> 3.2	38.4	-0.359495901	0.191964723
> 6.4	38.4	0.548288806	-1.011825356
> 9.6	38.4	-0.597741695	0.046560408
> 12.8	38.4	-1.003803022	0.319904274
> 16	38.4	0.060157729	0.824035747
> 19.2	38.4	0.251793973	0.632946256
> 22.4	38.4	-1.148625922	1.257360463
> 25.6	38.4	-0.957779058	0.30201831
> 3.2	41.6	0.267364885	1.139304073
> 6.4	41.6	0.171656399	0.623288738
> 9.6	41.6	-0.912461925	0.820049388
> 12.8	41.6	-0.645228334	1.041375034
> 16	41.6	0.83148026	2.914000894
> 19.2	41.6	-1.187664129	2.756035789
> 22.4	41.6	-0.467591981	0.467880919
> 25.6	41.6	-0.121964677	0.878460826
> 3.2	44.8	-0.158281328	-0.266451488
> 6.4	44.8	-0.235932095	1.440612621
> 12.8	44.8	-0.436274739	-0.298237058
> 16	44.8	-1.136811015	-0.688794913
> 19.2	44.8	0.87349541	0.449964869
> 22.4	44.8	-0.661067533	-0.08882536
> 25.6	44.8	0.301101862	-0.255259598
> 3.2	48	1.265698729	0.948786211
> 6.4	48	-0.806040053	0.352201302
> 9.6	48	1.006815543	1.12999254
> 12.8	48	-0.0132082	0.581935908
> 16	48	-1.214201235	1.431842632
> 19.2	48	-0.323473333	-0.083379843
> 22.4	48	-0.003250658	-0.248324838
> 25.6	48	-0.168398809	-0.916859761
> 3.2	51.2	-0.854793312	-0.231928117
> 6.4	51.2	0.507014021	-0.49107152
> 9.6	51.2	3.104090518	-2.654897126
> 12.8	51.2	-0.193785381	0.067244345
> 16	51.2	-0.71632415	-1.46047875
> 19.2	51.2	2.480632323	-0.336852087
> 22.4	51.2	-0.534663506	-1.093553232
> 25.6	51.2	-0.124544351	-0.603005465
> 
> Using the following code I compute a negative cross-variogram:
> 
> dtrd = as.data.frame(dtrd)
> # create a grid onto which i will interpolate # first get a range in data x.range <- as.integer(range(dtrd[,1])) y.range <- as.integer(range(dtrd[,2])) # now expand to a grid with desired spacing grid <- expand.grid(x=seq(from=x.range[1], to=x.range[2], by=0.4), y=seq(from=y.range[1], to=y.range[2], by=1.6))
> plot(grid)
> 
> # Convert to spatial data frame
> coordinates(grid) = c("x", "y")
> coordinates(dtrd) = c("x", "y")
> 
> # Compute cross-variogram
> g <- gstat(id = "n.ae", formula = n.ae ~ + 1, data = dtrd, nmax=10) g <- gstat(g, "n.ts", n.ts ~ 1, dtrd, nmax=10)
> g1 <- gstat(g, id=c("n.ae", "n.ts"), model = vgm(psill=-0.65, model="Sph", range=15, nugget=-0.05), fill.all = T) v <- variogram(g) g.fit<-fit.lmc(v, g1, correct.diagonal = 1.01)
> 
> # Plot cross-variogram with fitted model plot(v[v$id == "n.ae.n.ts",], g.fit$model[[2]], main="Cross-variogram: 
> Air Entry and Saturated Water Content", xlab="Lag Distance (mm)", ylab="Semivariance (-)", ylim=c(0,-0.8))
> 
> When attempting to use the negative cross-variogram in the co-kriging plan I receive the following error message:
> 
> # Perform co-kriging predictions
> ck <- predict.gstat(g1, grid)
> 
> non-positive definite coefficient matrix in structure 1non-positive definite coefficient matrix in structure 2Warning: No Intrinsic Correlation or Linear Model of Coregionalization found
> Reason: coefficient matrix not positive definite
> Warning: [add `set nocheck = 1;' to the command file to ignore the following error]
> 
> Error in predict.gstat(g1, grid) : 
>   gstat: value not allowed for: variograms do not satisfy a legal model
> 
> I have successfully produced kriged maps of each of the parameters individually, but do not know exactly how to deal with the negative cross-variogram and resolve the issue of a non positive definite correlation matrix. I think the problem arises from using a negative sill and nugget to model the sample cross-variogram with the fit.lmc function. If I use the same code but make the sill and nugget positive I get the same fit of the model to the sample cross-variogram and a message reading : Intrinsic Correlation found. Good. [using ordinary cokriging], which I initial thought that this was a good sign, but once again I can not resolve the following error message:
> 
> # Compute cross-variogram
> g <- gstat(id = "n.ae", formula = n.ae ~ + 1, data = dtrd, nmax=10) g <- gstat(g, "n.ts", n.ts ~ 1, dtrd, nmax=10)
> g1 <- gstat(g, id=c("n.ae", "n.ts"), model = vgm(psill=0.65, model="Sph", range=15, nugget=0.05), fill.all = T) v <- variogram(g) g.fit<-fit.lmc(v, g1, correct.diagonal = 1.01)
> 
> # Perform co-kriging predictions
> ck <- predict.gstat(g1, grid)
> 
> Intrinsic Correlation found. Good. 
> [using ordinary cokriging]
> 
> ERROR: "chfactor.c", line 131: singular matrix in function LDLfactor() Error in predict.gstat(g1, grid) : LDLfactor
> 
> The gstat manual shows that this problem can be fixed by adding correct.diagonal = 1.01 to the fit.lmc function, but this does not eliminate the error message
> 
> I am very confused about all of this and would appreciate any help to fix the non positive correlation matrix so that the negative cross-variogram can be used to perform co-kriging.
> 
> Thanks in advance for any help,
> 
> Blake Roberts
> 
> _______________________________________________
> R-sig-Geo mailing list
> R-sig-Geo at r-project.org
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
> 

-- 
Edzer Pebesma,  Co-Editor-in-Chief Computers & Geosciences
Institute for Geoinformatics (ifgi), University of Münster
Heisenbergstraße 2, 48149 Münster, Germany. Phone: +49 251
83 33081 http://ifgi.uni-muenster.de GPG key ID 0xAC227795

-------------- next part --------------
A non-text attachment was scrubbed...
Name: signature.asc
Type: application/pgp-signature
Size: 490 bytes
Desc: OpenPGP digital signature
URL: <https://stat.ethz.ch/pipermail/r-sig-geo/attachments/20150107/ce8ba836/attachment.bin>