[GRASS-stats] another gstat error

Sat Jun 20 07:32:10 EDT 2009


Roger Bivand pisze:
> On Sat, 20 Jun 2009, Jarek Jasiewicz wrote:
>
>>
>>
>> Roger Bivand pisze:
>>> On Fri, 19 Jun 2009, Jarosław Jasiewicz wrote:
>>>
>>>> Hi
>>>> Another gstat error:
>>>>
>>>> var.krigged=krige(z~x+y,gr,grid,model=var.fitted)
>>>> [using universal kriging]
>>>>
>>>> "chfactor.c", line 130: singular matrix in function LDLfactor()
>>>> Error in predict.gstat(g, newdata = newdata, block = block, nsim = 
>>>> nsim, :
>>>> LDLfactor
>>>
>>> Jarek:
>>>
>>> Shall I repeat what I wrote this morning? 
>> ????
>
> My reply is in the list archives and nabble, so I don't know why you 
> didn't get it.
yes, I found it on archive, I don't know why I did'n recive it, sorry.
>
>>> Did you actually try omitting the x+y in the formula? 
>> Yes, before I send post I tried firrst z~(x+y)^2 and if it failed I 
>> used z~1 to try without drift. The same error ocurs
>
> Googling 'gstat "chfactor.c", line 130: singular matrix in function 
> LDLfactor', I got to:
>
> http://osdir.com/ml/lang.r.geo/2006-07/msg00049.html
>
> which says:
>
> HINT: Read the manual at http://www.gstat.org/ ;
> look for: Trouble shooting -> Error messages -> From meschach
>
> "chfactor.c", line 130: singular matrix in function LDLfactor()
> Error in predict.gstat(g, newdata = newdata, block = block, nsim = 
> nsim, :
> matrix library error: gstat: matrix library error: singular matrix
> ---end R output
>
> after consulting the gstat manual I found the following reference:
>
> "These two error messages may occur (a) during simulation, when an
> observation falls almost exactly at a simulation location (b) when two
> observations occur at identical location occur and noaverage was 
> defined in
> its data definition. Solution to (a): increase the value of zero, to (b):
> remove the noaverage. Read also the next section."
>
> My input file was converted from lat/long to UTM using spTransform. The
> sampling density is quite high and it is possible that after conversion
> from latlong to UTM because of rounding some measurements fall on 
> identical
> coordinates. Is this a possible explanation for the error?
>
> The GSTAT online manual suggests to remove "noaverage" or increase "zero"
> to overcome the problem. However, I have not figured out how to do this
> when running gstat under R.
>
> So it looks as though you may have points that are seen as identical 
> in your data set.
>
> Hope this helps,
>
> Roger
Yes, it is possible, I must find method to identyfy and remove 
duplicates I will try SPATSTAT's remove duplicate functions

Thanks an many
Jarek
>
>>> Including a drift in untransformed coordinates is very likely to 
>>> lead to numerical chaos. If z~1 works, then the smoking gun is the 
>>> x+y, isn't it? If you want them, transform them (offset with for 
>>> example elide() in maptools and transform back afterwards).
>>>
>>> Roger
>>>
>>>>
>>>> I can't imagine the reason
>>>>
>>>> regards
>>>> Jarek
>>>>
>>>> here is structure of my gr spatialPointsDataFrame used in analysis:
>>>>
>>>> str(gr)
>>>> Formal class 'SpatialPointsDataFrame' [package "sp"] with 5 slots
>>>> ..@ data       :'data.frame':    2201 obs. of  5 variables:
>>>> .. ..$ cat: num [1:2201] 2 3 4 5 6 7 8 9 10 11 ...
>>>> .. ..$ naz: Factor w/ 874 levels "1000","1001",..: 874 1 2 3 4 5 6 
>>>> 7 8 9 ...
>>>> .. ..$ x  : num [1:2201] 360175 360176 360178 360169 360167 ...
>>>> .. ..$ y  : num [1:2201] 512736 512736 512733 512719 512712 ...
>>>> .. ..$ z  : num [1:2201] 86.5 87.4 87.7 86.9 86 ...
>>>> ..@ coords.nrs : num(0)
>>>> ..@ coords     : num [1:2201, 1:2] 360175 360176 360178 360169 
>>>> 360167 ...
>>>> .. ..- attr(*, "dimnames")=List of 2
>>>> .. .. ..$ : NULL
>>>> .. .. ..$ : chr [1:2] "coords.x1" "coords.x2"
>>>> ..@ bbox       : num [1:2, 1:2] 359907 512481 360492 512778
>>>> .. ..- attr(*, "dimnames")=List of 2
>>>> .. .. ..$ : chr [1:2] "coords.x1" "coords.x2"
>>>> .. .. ..$ : chr [1:2] "min" "max"
>>>> ..@ proj4string:Formal class 'CRS' [package "sp"] with 1 slots
>>>> .. .. ..@ projargs: chr " +proj=tmerc +lat_0=0 +lon_0=19 +k=0.9993 
>>>> +x_0=500000 +y_0=-5300000 +no_defs +a=6378137 +rf=298.257222101 
>>>> +towgs84=0.000,0.000,"| __truncated__
>>>>
>>>> and newdata grid:
>>>> str(grid)
>>>> Formal class 'SpatialGridDataFrame' [package "sp"] with 6 slots
>>>> ..@ data       :'data.frame':    131352 obs. of  1 variable:
>>>> .. ..$ region2: int [1:131352] NA NA NA NA NA NA NA NA NA NA ...
>>>> ..@ grid       :Formal class 'GridTopology' [package "sp"] with 3 
>>>> slots
>>>> .. .. ..@ cellcentre.offset: Named num [1:2] 359902 512478
>>>> .. .. .. ..- attr(*, "names")= chr [1:2] "x" "y"
>>>> .. .. ..@ cellsize         : num [1:2] 1.001 0.998
>>>> .. .. ..@ cells.dim        : int [1:2] 421 312
>>>> ..@ grid.index : int(0)
>>>> ..@ coords     : num [1:2, 1:2] 359902 360322 512478 512788
>>>> .. ..- attr(*, "dimnames")=List of 2
>>>> .. .. ..$ : NULL
>>>> .. .. ..$ : chr [1:2] "x" "y"
>>>> ..@ bbox       : num [1:2, 1:2] 359901 512477 360323 512789
>>>> .. ..- attr(*, "dimnames")=List of 2
>>>> .. .. ..$ : chr [1:2] "x" "y"
>>>> .. .. ..$ : chr [1:2] "min" "max"
>>>> ..@ proj4string:Formal class 'CRS' [package "sp"] with 1 slots
>>>> .. .. ..@ projargs: chr " +proj=tmerc +lat_0=0 +lon_0=19 +k=0.9993 
>>>> +x_0=500000 +y_0=-5300000 +no_defs +a=6378137 +rf=298.257222101 
>>>> +towgs84=0.000,0.000,"| __truncated__
>>>> _______________________________________________
>>>> grass-stats mailing list
>>>> grass-stats at lists.osgeo.org
>>>> http://lists.osgeo.org/mailman/listinfo/grass-stats
>>>>
>>>
>>
>