Re: [R-sig-Geo] Assessing and ranking the relationships and contribution of environmental correlates to species richness
So, it seems you have one raster for your dependent variable: species
richness and one for each explanatory variable.
I have never accomplished this but I think you should make a raster stack of
all these variables using raster package: the dependent one + the
explanatory ones then fit a model, for instance using a linear regression
(lm function)
With the model, you'll get coefficients for each exp. variables and
p-values. You'll trhen see which ones are significant and contribute the
most to explain species richness.
Hope this helps..
Mathieu
2011/7/4 ah3881
> Hi
>
> Thanks.
> I have species richness projections (made using Maxent projections of
> species distribution of 171 species, reclassified to give binary estimates
> of species distributions {10 percentile training threshold) before
> combining
> all 171 projections).
> So I have a raster of species richness in addition to all the possible
> correlates I am examining (npp, intra-annual variability of npp,
> inter-annual variability of npp, distance from coast, distance from karsts,
> latitude, lgm minimum temperature, lgm mean temperature, temperature change
> in minimum temperature since lgm, t change in mean temp since lgm, actual
> evapotranspiration)
> None of the correlates were used in species distribution projections and I
> have a gis raster of each variable, in addition to a database denoting the
> value of each for every km2 throughout Southeast Asia.
>
> Now I want to determine how much variability in species richness can be
> explained by variation in each of the other variables, in addition to the
> relationships. But I don't know which test would be best to do this, or
> software (GIS or database) which could best analyse this-especially given
> that the database is over 452 rows, and a large number of columns.
> Advice would be greatly appreciated
>
> --
> View this message in context:
> http://r-sig-geo.2731867.n2.nabble.com/Assessing-and-ranking-the-relationships-and-contribution-of-environmental-correlates-to-species-richs-tp6543945p6545267.html
> Sent from the R-sig-geo mailing list archive at Nabble.com.
>
> ___
> R-sig-Geo mailing list
> R-sig-Geo@r-project.org
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>
[[alternative HTML version deleted]]
___
R-sig-Geo mailing list
R-sig-Geo@r-project.org
https://stat.ethz.ch/mailman/listinfo/r-sig-geo
[R-sig-Geo] raster: cannot extract object from list
A function of mine outputs several raster objects in a list "l":
> class(l$lisaobsv)
[1] "RasterLayer"
attr(,"package")
[1] "raster"
But by some reason I get this weird error:
> a = l$lisaobsv
Error in slot(x, slotname) :
no slot of name "z" for this object of class "RasterLayer"
Calls: -> lapply -> FUN -> slot -> .Call
Warning message:
In .rk.get.structure.global("a") : failure to get object a
Is this a bug or is there an special way to extract raster objects from a list?
Thanks
Agus
> sessionInfo()
R version 2.13.0 (2011-04-13)
Platform: x86_64-pc-linux-gnu (64-bit)
locale:
[1] LC_CTYPE=en_US.UTF-8 LC_NUMERIC=C
[3] LC_TIME=en_US.UTF-8 LC_COLLATE=en_US.UTF-8
[5] LC_MONETARY=en_US.UTF-8 LC_MESSAGES=en_US.UTF-8
[7] LC_PAPER=en_US.UTF-8 LC_NAME=en_US.UTF-8
[9] LC_ADDRESS=en_US.UTF-8LC_TELEPHONE=en_US.UTF-8
[11] LC_MEASUREMENT=en_US.UTF-8LC_IDENTIFICATION=en_US.UTF-8
attached base packages:
[1] grid stats graphics grDevices utils datasets methods
[8] base
other attached packages:
[1] rgdal_0.6-33 geoR_1.6-27 plotrix_3.2-2 ggplot2_0.8.9 proto_0.3-8
[6] reshape_0.8.4 plyr_1.4 raster_1.8-40 sp_0.9-76 rkward_0.5.6
loaded via a namespace (and not attached):
[1] lattice_0.19-26 tools_2.13.0
___
R-sig-Geo mailing list
R-sig-Geo@r-project.org
https://stat.ethz.ch/mailman/listinfo/r-sig-geo
Re: [R-sig-Geo] raster: cannot extract object from list
Solved: it was caused by having upgraded to raster 1.8.40 but using
raster objects created under 1.8.38
Agus
2011/7/6 Agustin Lobo :
> A function of mine outputs several raster objects in a list "l":
>> class(l$lisaobsv)
> [1] "RasterLayer"
> attr(,"package")
> [1] "raster"
>
> But by some reason I get this weird error:
>> a = l$lisaobsv
> Error in slot(x, slotname) :
> no slot of name "z" for this object of class "RasterLayer"
> Calls: -> lapply -> FUN -> slot -> .Call
> Warning message:
> In .rk.get.structure.global("a") : failure to get object a
>
> Is this a bug or is there an special way to extract raster objects from a
> list?
>
> Thanks
>
> Agus
>
>> sessionInfo()
> R version 2.13.0 (2011-04-13)
> Platform: x86_64-pc-linux-gnu (64-bit)
>
> locale:
> [1] LC_CTYPE=en_US.UTF-8 LC_NUMERIC=C
> [3] LC_TIME=en_US.UTF-8 LC_COLLATE=en_US.UTF-8
> [5] LC_MONETARY=en_US.UTF-8 LC_MESSAGES=en_US.UTF-8
> [7] LC_PAPER=en_US.UTF-8 LC_NAME=en_US.UTF-8
> [9] LC_ADDRESS=en_US.UTF-8 LC_TELEPHONE=en_US.UTF-8
> [11] LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=en_US.UTF-8
>
> attached base packages:
> [1] grid stats graphics grDevices utils datasets methods
> [8] base
>
> other attached packages:
> [1] rgdal_0.6-33 geoR_1.6-27 plotrix_3.2-2 ggplot2_0.8.9 proto_0.3-8
> [6] reshape_0.8.4 plyr_1.4 raster_1.8-40 sp_0.9-76 rkward_0.5.6
>
> loaded via a namespace (and not attached):
> [1] lattice_0.19-26 tools_2.13.0
>
___
R-sig-Geo mailing list
R-sig-Geo@r-project.org
https://stat.ethz.ch/mailman/listinfo/r-sig-geo
Re: [R-sig-Geo] Subsetting RasterBrick very slow
Dear Robert, Thank you for sorting this out. And thanks to all the others for their valuable comments. Christian -- View this message in context: http://r-sig-geo.2731867.n2.nabble.com/Subsetting-RasterBrick-very-slow-tp6537890p6554009.html Sent from the R-sig-geo mailing list archive at Nabble.com. ___ R-sig-Geo mailing list R-sig-Geo@r-project.org https://stat.ethz.ch/mailman/listinfo/r-sig-geo
[R-sig-Geo] Summary: how to plot several raster layers on the same page and with the same color scale
I summarize here the different solutions I got to my question regarding
how to plot several raster layers on the same page and with the same
color scale.
Sorry it took a long time, I was interrupted by something else for weeks.
Hope this is useful.
Many thanks to all.
require(raster)
require(ggplot2)
require(plotrix)
require(rasterVis)
require(lattice)
foo1 <- foo2 <- foo3<- foo4 <- raster(nrows=20, ncols=20)
values(foo1) <- runif(400,-100.100)
values(foo2) <- runif(400,-150,100)
values(foo3) <- runif(400,0,50)
values(foo4) <- runif(400,-100,150)
extent(foo1) = extent(foo2) = extent(foo3) = extent(foo4) = extent(c(0,20,0,20))
#Based on comment by Jeremy Raw and andy.b...@wwu.edu
par(mfrow=c(2,2))
plot(foo1,col=rev(heat.colors(32)),zlim=c(-150,150))
plot(foo2,col=rev(heat.colors(32)),zlim=c(-150,150))
plot(foo3,col=rev(heat.colors(32)),zlim=c(-150,150))
plot(foo4,col=rev(heat.colors(32)),zlim=c(-150,150))
#Note by ALOBO: simplest solution
#Based on comment by Adam Sparks
require(plotrix)
minValue = -150
maxValue = 150
colstep=10
brk = c(minValue, maxValue, by = colstep)
brk = seq(minValue, maxValue, by = colstep)
par(mfrow=c(2,2))
plot(foo1, breaks = brk, col = rev(heat.colors(30)), axes = FALSE,
xlab = '', ylab = '', legend = FALSE)
color.legend(xl=20.5,yb=1,xr=21.5,yt=18, legend = seq(minValue,
maxValue, by = 3*colstep),
rect.col = rev(heat.colors(30)),pos=4,gradient="y",cex=0.55,offset=2,
col="black")
plot(foo2, breaks = brk, col = rev(heat.colors(30)),legend = FALSE,
xlab = '', ylab = '',)
color.legend(xl=20.5,yb=1,xr=21.5,yt=18, legend = seq(minValue,
maxValue, by = 3*colstep),
rect.col = rev(heat.colors(30)),pos=4,gradient="y",cex=0.55,offset=2,
col="black")
plot(foo3, breaks = brk, col = rev(heat.colors(30)),legend = FALSE,
xlab = '', ylab = '',)
color.legend(xl=20.5,yb=1,xr=21.5,yt=18, legend = seq(minValue,
maxValue, by = 3*colstep),
rect.col = rev(heat.colors(30)),pos=4,gradient="y",cex=0.55,offset=2,
col="black")
plot(foo3, breaks = brk, col = rev(heat.colors(30)),legend = FALSE,
xlab = '', ylab = '',)
color.legend(xl=20.5,yb=1,xr=21.5,yt=18, legend = seq(minValue,
maxValue, by = 3*colstep),
rect.col = rev(heat.colors(30)),pos=4,gradient="y",cex=0.55,offset=2,
col="black")
#Note by ALOBO: in practice, you must use hist() to determine
meaningful min and max values in both cases
#Based on comment by paul.hiems...@knmi.nl
library(ggplot2)
library(raster)
theme_set(theme_bw())
r = foo1
# Convert to data.frame
r_df = as.data.frame(as(r, 'SpatialPixelsDataFrame'))
print(str(r_df))
# create new column 10 times larger
# here you could extract data from other grids that
# you want to show at the same time
r_df$foo2 = values(foo2)
r_df$foo3 = values(foo3)
r_df$foo4 = values(foo4)
print(str(r_df))
# Reshape the data for ggplot
plotData = melt(r_df, id.vars = c('x','y'))
ggplot(aes(x = x, y = y), data = plotData) +
geom_tile(aes(fill = value)) + facet_wrap(~ variable) +
scale_fill_gradient(low = 'white', high = 'blue') +
coord_equal()
#Note by ALOBO: while this is the most aesthetic solution, would it
work for large raster objects?
#Based on comment by Robert Hijmans:
#Requires raster version >= 1.8-33 and rasterVis
s <- stack(foo1,foo2, foo3,foo4)
spplot(s)
theme_set(theme_bw())
gplot(s) + geom_tile(aes(fill = value)) + facet_wrap(~ variable) +
scale_fill_gradient(low = 'white', high = 'blue') + coord_equal()
#Note by ALOBO: assumed it automatically makes a sampleRegular() for
large raster objects?
#In such a case, it would combine aesthetic and performance
___
R-sig-Geo mailing list
R-sig-Geo@r-project.org
https://stat.ethz.ch/mailman/listinfo/r-sig-geo
Re: [R-sig-Geo] Summary: how to plot several raster layers on the same page and with the same color scale
Hi,
Just to add one more solution. It uses the rasterVis package (at
R-Forge now and should be at CRAN today).
library(raster)
library(rasterVis)
levelplot(s, par.settings=RdBuTheme)
More examples here: http://rastervis.r-forge.r-project.org/
Best,
Oscar.
---
Oscar Perpiñán Lamigueiro
Dpto. Ingeniería Eléctrica
EUITI-UPM
http://procomun.wordpress.com
-
>I summarize here the different solutions I got to my question regarding
>how to plot several raster layers on the same page and with the same
>color scale.
>Sorry it took a long time, I was interrupted by something else for
>weeks. Hope this is useful.
>Many thanks to all.
>
>
>require(raster)
>require(ggplot2)
>require(plotrix)
>require(rasterVis)
>require(lattice)
>
>foo1 <- foo2 <- foo3<- foo4 <- raster(nrows=20, ncols=20)
>values(foo1) <- runif(400,-100.100)
>values(foo2) <- runif(400,-150,100)
>values(foo3) <- runif(400,0,50)
>values(foo4) <- runif(400,-100,150)
>extent(foo1) = extent(foo2) = extent(foo3) = extent(foo4) =
>extent(c(0,20,0,20))
>
>
>#Based on comment by Jeremy Raw and andy.b...@wwu.edu
>par(mfrow=c(2,2))
>plot(foo1,col=rev(heat.colors(32)),zlim=c(-150,150))
>plot(foo2,col=rev(heat.colors(32)),zlim=c(-150,150))
>plot(foo3,col=rev(heat.colors(32)),zlim=c(-150,150))
>plot(foo4,col=rev(heat.colors(32)),zlim=c(-150,150))
>
>#Note by ALOBO: simplest solution
>
>#Based on comment by Adam Sparks
>require(plotrix)
>minValue = -150
>maxValue = 150
>colstep=10
>brk = c(minValue, maxValue, by = colstep)
>brk = seq(minValue, maxValue, by = colstep)
>par(mfrow=c(2,2))
>plot(foo1, breaks = brk, col = rev(heat.colors(30)), axes = FALSE,
>xlab = '', ylab = '', legend = FALSE)
>color.legend(xl=20.5,yb=1,xr=21.5,yt=18, legend = seq(minValue,
>maxValue, by = 3*colstep),
> rect.col = rev(heat.colors(30)),pos=4,gradient="y",cex=0.55,offset=2,
>col="black")
>plot(foo2, breaks = brk, col = rev(heat.colors(30)),legend = FALSE,
>xlab = '', ylab = '',)
>color.legend(xl=20.5,yb=1,xr=21.5,yt=18, legend = seq(minValue,
>maxValue, by = 3*colstep),
> rect.col = rev(heat.colors(30)),pos=4,gradient="y",cex=0.55,offset=2,
>col="black")
>plot(foo3, breaks = brk, col = rev(heat.colors(30)),legend = FALSE,
>xlab = '', ylab = '',)
>color.legend(xl=20.5,yb=1,xr=21.5,yt=18, legend = seq(minValue,
>maxValue, by = 3*colstep),
> rect.col = rev(heat.colors(30)),pos=4,gradient="y",cex=0.55,offset=2,
>col="black")
>plot(foo3, breaks = brk, col = rev(heat.colors(30)),legend = FALSE,
>xlab = '', ylab = '',)
>color.legend(xl=20.5,yb=1,xr=21.5,yt=18, legend = seq(minValue,
>maxValue, by = 3*colstep),
> rect.col = rev(heat.colors(30)),pos=4,gradient="y",cex=0.55,offset=2,
>col="black")
>
>
>#Note by ALOBO: in practice, you must use hist() to determine
>meaningful min and max values in both cases
>
>
>#Based on comment by paul.hiems...@knmi.nl
>library(ggplot2)
>library(raster)
>theme_set(theme_bw())
>
>r = foo1
>
># Convert to data.frame
>r_df = as.data.frame(as(r, 'SpatialPixelsDataFrame'))
>print(str(r_df))
># create new column 10 times larger
># here you could extract data from other grids that
># you want to show at the same time
>r_df$foo2 = values(foo2)
>r_df$foo3 = values(foo3)
>r_df$foo4 = values(foo4)
>
>print(str(r_df))
>
># Reshape the data for ggplot
>plotData = melt(r_df, id.vars = c('x','y'))
>
>ggplot(aes(x = x, y = y), data = plotData) +
>geom_tile(aes(fill = value)) + facet_wrap(~ variable) +
>scale_fill_gradient(low = 'white', high = 'blue') +
>coord_equal()
>
>#Note by ALOBO: while this is the most aesthetic solution, would it
>work for large raster objects?
>
>
>#Based on comment by Robert Hijmans:
>#Requires raster version >= 1.8-33 and rasterVis
>
>s <- stack(foo1,foo2, foo3,foo4)
>spplot(s)
>theme_set(theme_bw())
>gplot(s) + geom_tile(aes(fill = value)) + facet_wrap(~ variable) +
> scale_fill_gradient(low = 'white', high = 'blue') +
> coord_equal()
>
>#Note by ALOBO: assumed it automatically makes a sampleRegular() for
>large raster objects?
>#In such a case, it would combine aesthetic and performance
___
R-sig-Geo mailing list
R-sig-Geo@r-project.org
https://stat.ethz.ch/mailman/listinfo/r-sig-geo
Re: [R-sig-Geo] plotting latlon axes of sp object of different projection
Horacio See this thread from December: http://www.mail-archive.com/r-sig-geo@r-project.org/msg00154.html I haven't made much progress or further tested this approach since then. It seems to work for me and my needs. The example with the North Carolina map relies on pretty() to determine where the lines are. Depending on your projection, you may find it better to specify the graticule lines. cheers Josh On Jul 4, 2011, at 10:12 PM, Horacio Samaniego wrote: > I need to add a lat/lon axes to an equal area plot of a > SpatialLinesDataFrame object currently in units of meters. Is that > achievable using sp/maptools? > > The .prj file where the object comes from is: > PROJCS["Lambert_Azimuthal_Equal_Area",GEOGCS["GCS_unnamed > ellipse",DATUM["D_unknown",SPHEROID["Unknown",6370997,0]],PRIMEM["Greenwich",0],UNIT["Degree",0.017453292519943295]],PROJECTION["Lambert_Azimuthal_Equal_Area"],PARAMETER["latitude_of_origin",45],PARAMETER["central_meridian",-100],PARAMETER["false_easting",0],PARAMETER["false_northing",0],UNIT["Meter",1]] > > I'ld hate to install arcinfo to generate 1 plot! (plus the plot is > already nice as it is!!! just need the axes!) > > thanks, > I'll summarize! > > H > > -- > Horacio Samaniego > > Center for Nonlinear Studies > Los Alamos National Laboratory > MS B258 > Los Alamos, New Mexico 87545 > > Instituto de Silvicultura > Universidad Austral de Chile > Valdivia, Chile > > ___ > R-sig-Geo mailing list > R-sig-Geo@r-project.org > https://stat.ethz.ch/mailman/listinfo/r-sig-geo Josh London Wildlife Biologist National Marine Mammal Laboratory NOAA Seattle, Washington ___ R-sig-Geo mailing list R-sig-Geo@r-project.org https://stat.ethz.ch/mailman/listinfo/r-sig-geo
[R-sig-Geo] sp::overlay() not conserving row.names ?
I'm using sp::overlay between a SpPolDF and a SpPointsDF: > delme <- overlay(rndp,q1km, fn = mean) > class(q1km) [1] "SpatialPolygonsDataFrame" attr(,"package") [1] "sp" > class(rndp) [1] "SpatialPointsDataFrame" attr(,"package") [1] "sp" > class(delme) [1] "data.frame" but I'm confused by the fact that the row.names(delme) are not the same as the row.names(q1km): > delme[1:10,] ID id val X1 356. 357. 643.3162 X2 454. 455. 653.1833 X3 305. 306. 101.8930 X4 486.7143 487.7143 448.4612 X5 456.2500 457.2500 544.1646 X6 255.4000 256.4000 139.7507 X7 403.6667 404.6667 659.8114 X8 462.5000 463.5000 537.7209 X9 629.2000 630.2000 399.1040 X10 416.5000 417.5000 278.2879 > q1km@data[1:10,] idname ID XMIN XMAX YMIN YMAX X0 X0 0 -8.994 -8.594 43.595 43.995 X1 X1 1 -8.594 -8.194 43.595 43.995 X2 X2 2 -8.194 -7.794 43.595 43.995 X3 X3 3 -7.794 -7.394 43.595 43.995 X4 X4 4 -7.394 -6.994 43.595 43.995 X5 X5 5 -6.994 -6.594 43.595 43.995 X6 X6 6 -6.594 -6.194 43.595 43.995 X7 X7 7 -6.194 -5.794 43.595 43.995 X8 X8 8 -5.794 -5.394 43.595 43.995 X9 X9 9 -5.394 -4.994 43.595 43.995 Is this the way it has to be? If it is, then it's kind of inconvenient as I would like to join the new table "delme" to q1km@data. It seems I just have to subtract 1 to the numeric part of row.names(delme), but would like to make sure. Agus ___ R-sig-Geo mailing list R-sig-Geo@r-project.org https://stat.ethz.ch/mailman/listinfo/r-sig-geo
Re: [R-sig-Geo] sp::overlay() not conserving row.names ?
Agus, sp::overlay will be deprecated at some stage in favour of sp::over, which is a more consistent and complete approach to the same problem. For instance, overlay(x,y) would do the same as overlay(y,x), which is not good, if you think about it a bit longer (my own mistake, long time ago). Could you please check if sp::over has the same problems? See also http://cran.r-project.org/web/packages/sp/vignettes/over.pdf On 07/06/2011 05:28 PM, Agustin Lobo wrote: > I'm using sp::overlay between a SpPolDF and a SpPointsDF: > >> delme <- overlay(rndp,q1km, fn = mean) >> class(q1km) > [1] "SpatialPolygonsDataFrame" > attr(,"package") > [1] "sp" >> class(rndp) > [1] "SpatialPointsDataFrame" > attr(,"package") > [1] "sp" >> class(delme) > [1] "data.frame" > > but I'm confused by the fact that the row.names(delme) are not the same > as the row.names(q1km): >> delme[1:10,] > ID id val > X1 356. 357. 643.3162 > X2 454. 455. 653.1833 > X3 305. 306. 101.8930 > X4 486.7143 487.7143 448.4612 > X5 456.2500 457.2500 544.1646 > X6 255.4000 256.4000 139.7507 > X7 403.6667 404.6667 659.8114 > X8 462.5000 463.5000 537.7209 > X9 629.2000 630.2000 399.1040 > X10 416.5000 417.5000 278.2879 >> q1km@data[1:10,] >idname ID XMIN XMAX YMIN YMAX > X0 X0 0 -8.994 -8.594 43.595 43.995 > X1 X1 1 -8.594 -8.194 43.595 43.995 > X2 X2 2 -8.194 -7.794 43.595 43.995 > X3 X3 3 -7.794 -7.394 43.595 43.995 > X4 X4 4 -7.394 -6.994 43.595 43.995 > X5 X5 5 -6.994 -6.594 43.595 43.995 > X6 X6 6 -6.594 -6.194 43.595 43.995 > X7 X7 7 -6.194 -5.794 43.595 43.995 > X8 X8 8 -5.794 -5.394 43.595 43.995 > X9 X9 9 -5.394 -4.994 43.595 43.995 > > Is this the way it has to be? If it is, then it's kind of inconvenient > as I would like to join the new table "delme" to q1km@data. > It seems I just have to subtract 1 to the numeric part of > row.names(delme), but would like to make sure. > > Agus > > ___ > R-sig-Geo mailing list > R-sig-Geo@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.pebe...@wwu.de ___ R-sig-Geo mailing list R-sig-Geo@r-project.org https://stat.ethz.ch/mailman/listinfo/r-sig-geo
Re: [R-sig-Geo] sp::overlay() not conserving row.names ?
Yes! I'm doing: > delme2 <- over(q1km, rndp,fn = mean) And now not only the row.names are correct and the order of the 2 arguments makes more sense, but the dimensions of the result are also the same that that of the first argument (which implies rows of NA for those polygons with no points: the user can get rid of these rows with na.omit() later, for me having the NA is a lot better): > delme2[1:3,] ID val 0 356. 760.5662 1 454. 541.5580 2 305. 467.0864 > q1km@data[1:3,] ID XMIN XMAX YMIN YMAX 0 0 -8.994 -8.594 43.595 43.995 1 1 -8.594 -8.194 43.595 43.995 2 2 -8.194 -7.794 43.595 43.995 > dim(delme2) [1] 299 3 > dim(q1km@data) [1] 299 6 I can just join this table to the one of q1km: > delme2 = cbind(idname=row.names(delme2),delme2) > delme2[1:3,] idname ID val 0 0 356. 760.5662 1 1 454. 541.5580 2 2 305. 467.0864 > q1kmv2 = q1km > q1kmv2 <- mijoin(q1kmv2,delme2[,-2], by.x=1,by.y=1) > q1kmv2@data[1:3,] ID XMIN XMAX YMIN YMAX idname val 0 0 -8.994 -8.594 43.595 43.995 0 760.5662 1 1 -8.594 -8.194 43.595 43.995 1 541.5580 2 2 -8.194 -7.794 43.595 43.995 2 467.0864 3 3 -7.794 -7.394 43.595 43.995 3 487.6461 Thanks! Agus 2011/7/6 Edzer Pebesma : > Agus, sp::overlay will be deprecated at some stage in favour of > sp::over, which is a more consistent and complete approach to the same > problem. For instance, overlay(x,y) would do the same as overlay(y,x), > which is not good, if you think about it a bit longer (my own mistake, > long time ago). > > Could you please check if sp::over has the same problems? See also > http://cran.r-project.org/web/packages/sp/vignettes/over.pdf > > On 07/06/2011 05:28 PM, Agustin Lobo wrote: >> I'm using sp::overlay between a SpPolDF and a SpPointsDF: >> >>> delme <- overlay(rndp,q1km, fn = mean) >>> class(q1km) >> [1] "SpatialPolygonsDataFrame" >> attr(,"package") >> [1] "sp" >>> class(rndp) >> [1] "SpatialPointsDataFrame" >> attr(,"package") >> [1] "sp" >>> class(delme) >> [1] "data.frame" >> >> but I'm confused by the fact that the row.names(delme) are not the same >> as the row.names(q1km): >>> delme[1:10,] >> ID id val >> X1 356. 357. 643.3162 >> X2 454. 455. 653.1833 >> X3 305. 306. 101.8930 >> X4 486.7143 487.7143 448.4612 >> X5 456.2500 457.2500 544.1646 >> X6 255.4000 256.4000 139.7507 >> X7 403.6667 404.6667 659.8114 >> X8 462.5000 463.5000 537.7209 >> X9 629.2000 630.2000 399.1040 >> X10 416.5000 417.5000 278.2879 >>> q1km@data[1:10,] >> idname ID XMIN XMAX YMIN YMAX >> X0 X0 0 -8.994 -8.594 43.595 43.995 >> X1 X1 1 -8.594 -8.194 43.595 43.995 >> X2 X2 2 -8.194 -7.794 43.595 43.995 >> X3 X3 3 -7.794 -7.394 43.595 43.995 >> X4 X4 4 -7.394 -6.994 43.595 43.995 >> X5 X5 5 -6.994 -6.594 43.595 43.995 >> X6 X6 6 -6.594 -6.194 43.595 43.995 >> X7 X7 7 -6.194 -5.794 43.595 43.995 >> X8 X8 8 -5.794 -5.394 43.595 43.995 >> X9 X9 9 -5.394 -4.994 43.595 43.995 >> >> Is this the way it has to be? If it is, then it's kind of inconvenient >> as I would like to join the new table "delme" to q1km@data. >> It seems I just have to subtract 1 to the numeric part of >> row.names(delme), but would like to make sure. >> >> Agus >> >> ___ >> R-sig-Geo mailing list >> R-sig-Geo@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.pebe...@wwu.de > > ___ > R-sig-Geo mailing list > R-sig-Geo@r-project.org > https://stat.ethz.ch/mailman/listinfo/r-sig-geo > ___ R-sig-Geo mailing list R-sig-Geo@r-project.org https://stat.ethz.ch/mailman/listinfo/r-sig-geo
Re: [R-sig-Geo] sp::overlay() not conserving row.names ?
I'm not too familiar with GIS terminology, but would you agree that sp::over does (or allows for) a spatial table join? On 07/06/2011 06:10 PM, Agustin Lobo wrote: > Yes! > > I'm doing: >> delme2 <- over(q1km, rndp,fn = mean) > > And now not only the row.names are correct and the order of the 2 > arguments makes more sense, but the dimensions of the result are also > the same > that that of the first argument > (which implies rows of NA for those polygons with no points: the user > can get rid of > these rows with na.omit() later, for me having the NA is a lot better): > >> delme2[1:3,] > ID val > 0 356. 760.5662 > 1 454. 541.5580 > 2 305. 467.0864 >> q1km@data[1:3,] > ID XMIN XMAX YMIN YMAX > 0 0 -8.994 -8.594 43.595 43.995 > 1 1 -8.594 -8.194 43.595 43.995 > 2 2 -8.194 -7.794 43.595 43.995 > >> dim(delme2) > [1] 299 3 >> dim(q1km@data) > [1] 299 6 > > I can just join this table to the one of q1km: >> delme2 = cbind(idname=row.names(delme2),delme2) >> delme2[1:3,] > idname ID val > 0 0 356. 760.5662 > 1 1 454. 541.5580 > 2 2 305. 467.0864 > >> q1kmv2 = q1km >> q1kmv2 <- mijoin(q1kmv2,delme2[,-2], by.x=1,by.y=1) >> q1kmv2@data[1:3,] > ID XMIN XMAX YMIN YMAX idname val > 0 0 -8.994 -8.594 43.595 43.995 0 760.5662 > 1 1 -8.594 -8.194 43.595 43.995 1 541.5580 > 2 2 -8.194 -7.794 43.595 43.995 2 467.0864 > 3 3 -7.794 -7.394 43.595 43.995 3 487.6461 > > Thanks! > Agus > > 2011/7/6 Edzer Pebesma : >> Agus, sp::overlay will be deprecated at some stage in favour of >> sp::over, which is a more consistent and complete approach to the same >> problem. For instance, overlay(x,y) would do the same as overlay(y,x), >> which is not good, if you think about it a bit longer (my own mistake, >> long time ago). >> >> Could you please check if sp::over has the same problems? See also >> http://cran.r-project.org/web/packages/sp/vignettes/over.pdf >> >> On 07/06/2011 05:28 PM, Agustin Lobo wrote: >>> I'm using sp::overlay between a SpPolDF and a SpPointsDF: >>> delme <- overlay(rndp,q1km, fn = mean) class(q1km) >>> [1] "SpatialPolygonsDataFrame" >>> attr(,"package") >>> [1] "sp" class(rndp) >>> [1] "SpatialPointsDataFrame" >>> attr(,"package") >>> [1] "sp" class(delme) >>> [1] "data.frame" >>> >>> but I'm confused by the fact that the row.names(delme) are not the same >>> as the row.names(q1km): delme[1:10,] >>> ID id val >>> X1 356. 357. 643.3162 >>> X2 454. 455. 653.1833 >>> X3 305. 306. 101.8930 >>> X4 486.7143 487.7143 448.4612 >>> X5 456.2500 457.2500 544.1646 >>> X6 255.4000 256.4000 139.7507 >>> X7 403.6667 404.6667 659.8114 >>> X8 462.5000 463.5000 537.7209 >>> X9 629.2000 630.2000 399.1040 >>> X10 416.5000 417.5000 278.2879 q1km@data[1:10,] >>>idname ID XMIN XMAX YMIN YMAX >>> X0 X0 0 -8.994 -8.594 43.595 43.995 >>> X1 X1 1 -8.594 -8.194 43.595 43.995 >>> X2 X2 2 -8.194 -7.794 43.595 43.995 >>> X3 X3 3 -7.794 -7.394 43.595 43.995 >>> X4 X4 4 -7.394 -6.994 43.595 43.995 >>> X5 X5 5 -6.994 -6.594 43.595 43.995 >>> X6 X6 6 -6.594 -6.194 43.595 43.995 >>> X7 X7 7 -6.194 -5.794 43.595 43.995 >>> X8 X8 8 -5.794 -5.394 43.595 43.995 >>> X9 X9 9 -5.394 -4.994 43.595 43.995 >>> >>> Is this the way it has to be? If it is, then it's kind of inconvenient >>> as I would like to join the new table "delme" to q1km@data. >>> It seems I just have to subtract 1 to the numeric part of >>> row.names(delme), but would like to make sure. >>> >>> Agus >>> >>> ___ >>> R-sig-Geo mailing list >>> R-sig-Geo@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.pebe...@wwu.de >> >> ___ >> R-sig-Geo mailing list >> R-sig-Geo@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.pebe...@wwu.de ___ R-sig-Geo mailing list R-sig-Geo@r-project.org https://stat.ethz.ch/mailman/listinfo/r-sig-geo
