[R-sig-Geo] spatialgridataframes

Lyndon Estes lestes at princeton.edu
Mon Feb 21 15:28:37 CET 2011


Hi Edzer,

Thanks for the clarification regarding factors.

Cheers, Lyndon

On Sun, Feb 20, 2011 at 3:50 PM, Edzer Pebesma
<edzer.pebesma at uni-muenster.de> wrote:
> Lyndon, the original question was about removing a particular level in a
> factor variable. The answer you gave seems to work for numeric variables
> only (funny enough, soil in meuse.grid is a numeric variable!), and
> package raster doesn't seem to deal well with factors (although coercion
> to raster does not complain). Building on your example:
>
> meuse.grid$soilf = factor(meuse.grid$soil, labels=c("cl1","cl2","cl3"))
> m2 <- as(meuse.grid, "SpatialGridDataFrame")
> m3 <- raster(m2, layer = "soilf")  # Convert soilf factor to raster
> plot(m3)
>
> It seems that here, the factor levels are lost, as well as the knowledge
> that this was a factor, at least it is lost when converted back to a
> SpatialGridDataFrame.
>
> Function factor can be used to change factors, or factor levels, or
> remove some, as in:
>
> meuse.grid$soilf2 = factor(meuse.grid$soilf,
> levels=levels(meuse.grid$soilf)[1:2])
>
> Working with factors converted to numeric values in linear regression
> models typically leads to plain wrong results.
>
> On 02/20/2011 05:44 PM, Lyndon Estes wrote:
>> Hi Mary,
>>
>> This might help answer your questions. I used the meuse dataset and
>> converted to raster formats, but I think the general approach should
>> work for what you want to do.
>>
>> library(raster)
>> library(gstat)
>> data(meuse.grid)
>> coordinates(meuse.grid) = ~x+y
>> gridded(meuse.grid) = TRUE
>> class(meuse.grid)
>> m2 <- as(meuse.grid, "SpatialGridDataFrame")
>> m3 <- raster(m2, layer = "soil")  # Convert soil classes to raster
>>
>>>
>>> Question1.
>>> How can i remove the ninth class in R because it does not have to be included in
>>> geostatical analysis.
>>
>> m4 <- m3 * ((m3 < 3) / (m3 < 3))  # Removes class 3 from soil,
>> converts it to NA values (this could also
>> # serve as a mask)
>>
>> # If you want to keep that part of the grid in the analysis, then you
>> might want to collapse the one class
>> # into another
>> m5 <- (m3 == 3 | m3 == 2) * 2 + (m3 == 1)  # Class 2 now includes 2 and 3
>>
>>> Question 2
>>> how can i use the map created in Q1 to clip the other 5  aboventioned maps (eg
>>> DEM etc) or how can i create a mask from the map in Q1.
>>>
>>
>>
>> # Create a mask for just the area of soil class 1
>> sc1.mask <- (m3 == 1) / (m3 == 1)
>>
>> # You then multiply your other rasters by your mask to reduce rasters
>> to the areas you want to analyze.
>>
>> Cheers, Lyndon
>>
>> _______________________________________________
>> R-sig-Geo mailing list
>> R-sig-Geo at r-project.org
>> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>
> --
> Edzer Pebesma
> Institute for Geoinformatics (ifgi), University of Münster
> Weseler Straße 253, 48151 Münster, Germany. Phone: +49 251
> 8333081, Fax: +49 251 8339763  http://ifgi.uni-muenster.de
> http://www.52north.org/geostatistics      e.pebesma at wwu.de
>
> _______________________________________________
> R-sig-Geo mailing list
> R-sig-Geo at r-project.org
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>



-- 
Lyndon Estes
Research Associate
Woodrow Wilson School
Princeton University
+1-609-258-2392 (o)
+1-609-258-6082 (f)
+1-202-431-0496 (m)
lestes at princeton.edu



More information about the R-sig-Geo mailing list