[GRASS-stats] Raster maps from GRASS to R and back to GRASS?

Dylan Beaudette debeaudette at ucdavis.edu
Tue Mar 3 11:53:43 EST 2009


On Tuesday 03 March 2009, Nikos Alexandris wrote:
> Hi GRASS&R experts,
>
> After playing around with R's functions for PCA, I need to
>
> 1. load GRASS raster maps in R
> 2. Do something with the data
>    -involved workaround for NA's
> 3. Get back the NA's
> 4. Write back the transformed data in GRASS as new raster maps.
>
>
> So, this is rather a (general?) GRASS - R question: how is it done,
> after having calculated the principal components for a given set of
> GRASS raster maps in R, to create new GRASS raster maps with them? I
> mean how is it done to write the transformed values back to the right
> coordinates?
>
> Is there a clear example on the web using the writeRAST6()function?
>

Does this help:
http://casoilresource.lawr.ucdavis.edu/drupal/node/664

Cheers,

Dylan


> Below an example.
> Regards, Nikos
> -------------------------------------------------------------------------
>
>  # launch R from within GRASS
> R
>
>
>  # load spgrass6, and projection parameters
> library(spgrass6)
> G <- gmeta6()
>
>
>  # load data
> modis_2006_250m <- readRAST6 (c('MOD2006_239_250_sur_refl_b01',
> 'MOD2006_239_250_sur_refl_b02'))
>
>
>  # workaround NA's as explained by Dylan [1]
> modis_2006_250m.values <- which(!is.na(modis_2006_250m at data
> $MOD2006_239_250_sur_refl_b01) & !is.na(modis_2006_250m at data
> $MOD2006_239_250_sur_refl_b02))
>
>
>  # create object without NA's
> mod06_250m <- modis_2006_250m at data[modis_2006_250m.values, ]
>
>
>  # perform pca based on SVD (default with center=TRUE, scale=FALSE)
> prcomp.mod06_250m <- prcomp(mod06_250m)
>
>
>  # check eigenvalues
> summary(prcomp.mod06_250m)
>
> Importance of components:
>                            PC1      PC2
> Standard deviation     3183.53 330.0584
> Proportion of Variance    0.99   0.0106
> Cumulative Proportion     0.99   1.0000
>
>
>  # check eigenvectors (=rotation)
> prcomp.mod06_250m$rotation
>
>                                     PC1        PC2
> MOD2006_239_250_sur_refl_b01 -0.3640966 -0.9313612
> MOD2006_239_250_sur_refl_b02 -0.9313612  0.3640966
>
>
>  # write back to GRASS the transformed data as new RASTER maps?
> --> Any help here on how to get the values to the right coordinates? Is
> it enough to just use some <--
>
>
>
>
> =======================================================================
> Extra info: comparing the resulted PC's from both R and GRASS
> =======================================================================
>  # compare with i.pca
> i.pca input="MOD06_250_b01,MOD06_250_b02" output="TEST06_250m_b12"
> rescale=0,0
>
> Eigen values, (vectors), and [percent importance]:
>    PC1 3481240.24 (  0.36  0.93 ) [ 98.92% ]
>    PC2  37875.46 (  0.93 -0.36 ) [  1.08% ]
>
>
>    ### Note ###
> The signs of the eigenvectors exactly the opposite. Therefore we expect
> that the Principal Components deriving from R's "prcomp()" will be the
> opposite of the Components deriving from "i.pca"
>    ### #### ###
>
>
>  # get summary for PC's produced by "prcomp()"
> summary(prcomp.mod06_250m$x)
>       PC1                 PC2
>  Min.   :-10167.54   Min.   :-3919.99
>  1st Qu.: -3422.88   1st Qu.: -167.49
>  Median : -3119.27   Median :   45.32
>  Mean   : -3135.72   Mean   :   21.60
>  3rd Qu.: -2851.33   3rd Qu.:  239.13
>  Max.   :   -64.22   Max.   : 1611.15
>
>
>  # get variances for PC's produce by "prcomp()"
> var(prcomp.mod06_250m$x[, 1])
> [1] 302107.8
>
> var(prcomp.mod06_250m$x[, 2])
> [1] 108471.9
>
>
>  # compare with PC's from GRASS
>
> ### kept are only commonly reported parameters by GRASS' "r.univar" and
> R's "summary()" ###
>
>
>  # for GRASS' PC1
> r.univar TEST2006_250m_b12.1
>
> [...]
> n: 350879
> minimum: 64.2261
> maximum: 10166.7
> ...
> mean: 3139.65
> mean of absolute values: 3139.65
> standard deviation: 543.404
> variance: 295287
>
>
>  # for GRASS' PC2
> r.univar TEST2006_250m_b12.2
>
>  100%
> total null and non-null cells: 1023328
> total null cells: 672449
>
> Of the non-null cells:
> ----------------------
> n: 350879
> minimum: -1609.93
> maximum: 3922.17
>
> mean: -22.0414
> ...
> variance: 108444
>
> ######## Result is #########
> ### almost identical :-) ###
> ############################
>
> _______________________________________________
> grass-stats mailing list
> grass-stats at lists.osgeo.org
> http://lists.osgeo.org/mailman/listinfo/grass-stats



-- 
Dylan Beaudette
Soil Resource Laboratory
http://casoilresource.lawr.ucdavis.edu/
University of California at Davis
530.754.7341


More information about the grass-stats mailing list