[GRASSLIST:2265] Re: trend surfaces - can't get image of residuals

[Roger Bivand]

> Hello all,
>
> Things are starting to work on trend surfaces analysis.
> I've converted my raster to sites and used sites.get() in R, so it comes
> as a data.frame, and surf.ls() works well. using trmat.G() and
> rast.put(), I can get the trend surface back into GRASS.

You have the fitted values (and the residuals) for the site locations used
to fit the model. In order to get fitted values for a raster grid, you
need to provide predict(b2, newdata= ...) with the raster cell locations,
most conveniently with trmat.G(G, b2), which gives the raster cell
predictions. The output can be visualsed in R with

plot(G, trmat.G(G, b2))

provided that the total number of cells is not too large. Or you can
choose to move it to GRASS with rast.put() and visualise it there. A lot
of this is covered in the examples at the end of the trmat.G help page.
See also Neteler and Mitasova "Open Source GIS", p. 361, and

http://grass.itc.it/statsgrass/grass_geostats.html

also has useful ideas.

Roger

>
> The problem now is with residuals. I don't know how to put these values
> into a matrix, to get an image or get it into GRASS.
>
> if I use residuals(), I can see thet its output is the same as the $wz
> in surf.ls() output...
>
> here goes some of what I have:
>
>> bases<-sites.get(G,slist="isobases",all.sites=F)
> Warning message:
> Old style sites file
>> str(bases)
> `data.frame':   3509 obs. of  4 variables:
>  $ east : num  678148 672602 672608 680412 663768 ...
>  $ north: num  7785498 7785472 7785468 7785468 7785462 ...
>  $ id   : int  321 322 323 324 325 326 327 328 329 330 ...
>  $ num1 : num  680 680 680 640 720 640 720 720 720 720 ...
>  - attr(*, "nncols")= Named int  2 -1 1 0
>   ..- attr(*, "names")= chr  "dims" "cattype" "dbl.att" "str.att"
>> b2<-surf.ls(2,bases$east,bases$north,bases$num1)
>
>> str(b2)
> List of 11
>  $ x   : num [1:3509] 678148 672602 672608 680412 663768 ...
>  $ y   : num [1:3509] 7785498 7785472 7785468 7785468 7785462 ...
>  $ z   : num [1:3509] 680 680 680 640 720 640 720 720 720 720 ...
>  $ np  : num 2
>  $ f   : num [1:21054] 1 1 1 1 1 1 1 1 1 1 ...
>  $ r   : num [1:21] -59.237   4.521 -36.100 -22.345   0.472 ...
>  $ beta: num [1:6]  945.2 -251.6  142.0   23.1  -48.5 ...
>  $ wz  : num [1:3509] 241.8 166.0 165.8 211.4 -46.8 ...
>  $ rx  : num [1:2] 657858 681048
>  $ ry  : num [1:2] 7762302 7785498
>  $ call: language surf.ls(np = 2, x = bases$east, y = bases$north, z =
> bases$num1)
>  - attr(*, "class")= chr "trls"
>
>> b2r<-residuals(b2)
>> str(b2r)
>  num [1:3509] 241.8 166.0 165.8 211.4 -46.8 ...
>
> thanks in advance
> --
> +-------------------------------------------------+
>         Carlos Henrique Grohmann - Guano
>     Geologist - MSc Student at IGc-USP - Brazil
>        Linux User #89721  ICQ: 214752832
> +-------------------------------------------------+
>
>
>> Hello all
>>
>> I'm trying to do some trend surface analysis with GRASS data in R, but
>> I with some problems:
>>
>> 1 - when I get the raster from GRASS, surf.ls() doesn't accept it,
>> since its not in x,y,z format:
>>
>>> str(isobase)
>> List of 1
>>  $ isobases.auto: num [1:241266] 771 770 769 768 767 ...
>>>
>>
>> So, how can I change this to a format that surf.ls() understands?
>>
>
> isobase$isobases.auto is the z argument. You can generate the x by
> east() of the object returned by gmeta(), and y by north() of this
> object. I'm not sure how well surf.ls() scales to this number of points,
> though.
>
>> 2 - I'll try several degrees of surfaces, so I'll use surf.ls(), but
>> what about residuals? can I get them with residuals()?
>>
> Yes, residuals() should work. Again, scaling to this number of points is
> an issue. On the other hand, the residuals and fitted values will drop
> straight into the GRASS location.
>
> Roger


-- 
Roger Bivand
NHH, Breiviksveien 40, N-5045 Bergen, Norway