Subscribe to R-sig-geo feed
This is an archive for R-sig-geo, not a forum. It is not possible to post through Nabble - you may not start a new thread nor follow up an existing thread. If you wish to post, but are not subscribed to the list through the list homepage, subscribe first (through the list homepage, not via Nabble) and post from your subscribed email address. Until 2015-06-20, subscribers could post through Nabble, but policy has been changed as too many non-subscribers misunderstood the interface.
Updated: 16 min 34 sec ago

Re: Spatial clustering with spdep::skater

Mon, 07/09/2018 - 07:22
Hi Axel,

The original paper suggests to look at the SSW decay. You can visualize
it easily as it is returned in the output from skater().

Best,

Elias


On 08/07/18 08:52, Axel Urbiz wrote:
> Dear Experts,
>
> My apologies in advance as this is more a statistical question, rather than an R programming question, but hope you can point me in the right direction.
>
> I’m working with spdep::skater to fit clusters to spatial data subject to contiguity constraints. This function fits clusters by edge removal from Minimum Spanning Trees.
>
> In this context, I’d appreciate any pointers on how to tune the number of clusters. What is a sensible criteria to use?
>
> Thanks,
> Axel.
> _______________________________________________
> R-sig-Geo mailing list
> [hidden email]
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo

Spatial clustering with spdep::skater

Sun, 07/08/2018 - 06:52
Dear Experts,

My apologies in advance as this is more a statistical question, rather than an R programming question, but hope you can point me in the right direction.

I’m working with spdep::skater to fit clusters to spatial data subject to contiguity constraints. This function fits clusters by edge removal from Minimum Spanning Trees.

In this context, I’d appreciate any pointers on how to tune the number of clusters. What is a sensible criteria to use?

Thanks,
Axel.
_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo

geom_contour

Fri, 07/06/2018 - 16:30
Hi, I'm trying to generate a contour plot from a data set I have in Excel
showing how a given set of three of the 8 columns/categories are related to
each other. I pulled out three arbitrary columns using the subset function,
but I'm kinda lost with what to do from there.

Lets say I have this subset:

      x       y       z
-------------------------
      7       6      5.3
      1     3.14   4.0
      5      2.7    6.1

How would I plot a contour graph of the subset? Is that even possible?
Thanks.

        [[alternative HTML version deleted]]

_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo

Re: Distance from not NA cells in a raster

Fri, 07/06/2018 - 12:47
Sorry, I forgot to include session info to associate my issue...thanks!
> sessionInfo()
R version 3.4.3 (2017-11-30)
Platform: x86_64-w64-mingw32/x64 (64-bit)
Running under: Windows 7 x64 (build 7601) Service Pack 1

other attached packages:
[1] RSAGA_1.2.0    plyr_1.8.4     shapefiles_0.7 foreign_0.8-69 gstat_1.1-6    rgdal_1.2-16
[7] raster_2.6-7   sp_1.2-5


From: Gregovich, Dave P (DFG)
Sent: Friday, July 06, 2018 8:23 AM
To: '[hidden email]'
Subject: Distance from not NA cells in a raster

Hi,
I would like to obtain the distance from all not NA cells in a raster. This works for smaller rasters, but seems difficult for the size of rasters (~ 8000 pixel square)  I am working with.
Below is what I've tried. I would be OK calling other software from R, or using some parallelization, if it might help.
Thanks so much for your help!  If I could just calculate this distance in two hours or less (or so) I would be satisfied.
Dave.

rm(list=ls())
library(raster)

#make raster
rast <- raster(nrow = 8000, ncol = 8000, ext = extent(0,1,0,1))

#generate cells to calculate distance from.
rast[sample(8000^2, 10000)] <- 1

#try two different methods...
dist1 <- gridDistance(rast, origin = 1)#throws an error after x minutes
#'Error: cannot allocate vector of size 3.8 Gb'
dist2 <- distance(rast)#ran all night, R was hung up in the morning and had to force shutdown.

___________________________________________
Dave Gregovich
Research Analyst
Alaska Department of Fish and Game
Division of Wildlife Conservation
802 3rd Street
Douglas, AK 99824
907-465-4291
___________________________________________


        [[alternative HTML version deleted]]

_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo

Distance from not NA cells in a raster

Fri, 07/06/2018 - 11:22
Hi,
I would like to obtain the distance from all not NA cells in a raster. This works for smaller rasters, but seems difficult for the size of rasters (~ 8000 pixel square)  I am working with.
Below is what I've tried. I would be OK calling other software from R, or using some parallelization, if it might help.
Thanks so much for your help!  If I could just calculate this distance in two hours or less (or so) I would be satisfied.
Dave.

rm(list=ls())
library(raster)

#make raster
rast <- raster(nrow = 8000, ncol = 8000, ext = extent(0,1,0,1))

#generate cells to calculate distance from.
rast[sample(8000^2, 10000)] <- 1

#try two different methods...
dist1 <- gridDistance(rast, origin = 1)#throws an error after x minutes
#'Error: cannot allocate vector of size 3.8 Gb'
dist2 <- distance(rast)#ran all night, R was hung up in the morning and had to force shutdown.

___________________________________________
Dave Gregovich
Research Analyst
Alaska Department of Fish and Game
Division of Wildlife Conservation
802 3rd Street
Douglas, AK 99824
907-465-4291
___________________________________________


        [[alternative HTML version deleted]]

_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo

Kriging and Spatio temporal kriging on air quality data

Thu, 07/05/2018 - 01:34
Hi, I am using kriging and spatio temporal kriging for studying and
analyzing air quality data from mobile sensors.
I would like to ask some question about the following issues:
1) is possible to gave the same values scale to different kriging result
(maps from different measures in the same study area). That is if I have a
kriging map with values in the range of 0-20 and another map with values in
the range of 5-30, can I gave to both the maps a values (and colours) scale
in the range of 0-30, in order to compare the values from the two maps.
2) is there some particular package and or tecnique specifically indicated
to make spatial kriging from air quality data from mobile sensor.

Kind regards.

        [[alternative HTML version deleted]]

_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo

Re: How to use gstat 'variogram' function for nonlinear formulae?

Wed, 07/04/2018 - 11:28


On 07/04/2018 06:14 PM, Joelle k. Akram wrote:
> Dear all,
>
>
> For a linear model I use the 'variogram' function in gstat as follows:
>
>
> coordinates(data) <- c('x','y')
> myvario <- variogram(dep~var1+var2+var3, data)
>
> Now, I have the following Nonlinear formula that I want to use in variogram.
>
> dep~ var1*var2^3
This is still a linear model.

>
>
> My 2 questions are:
>
> i) Can the variogram function in gstat handle nonlinear models as this formula?
>
> ii) What is the syntax to use this nonlinear formulae within the variogram function? is it simply :
>
> myvario <- variogram(dep~ var1*var2^3, data)

I think that like in other areas of R you'll have to use

myvario <- variogram(dep~ var1*I(var2)^3, data)

Note that the * will be interpreted as indicating you want main effects
and an interaction, not simply the product of var1*var3^3; see ?formula

>
> Any advice is appeciated.
>
>
> cheers
>
> Chris Akram
>
>
>
> [[alternative HTML version deleted]]
>
> _______________________________________________
> R-sig-Geo mailing list
> [hidden email]
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>
--
Edzer Pebesma
Institute for Geoinformatics
Heisenbergstrasse 2, 48151 Muenster, Germany
Phone: +49 251 8333081

_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo

How to use gstat 'variogram' function for nonlinear formulae?

Wed, 07/04/2018 - 11:14
Dear all,


For a linear model I use the 'variogram' function in gstat as follows:


coordinates(data) <- c('x','y')
myvario <- variogram(dep~var1+var2+var3, data)

Now, I have the following Nonlinear formula that I want to use in variogram.

dep~ var1*var2^3


My 2 questions are:

i) Can the variogram function in gstat handle nonlinear models as this formula?

ii) What is the syntax to use this nonlinear formulae within the variogram function? is it simply :

myvario <- variogram(dep~ var1*var2^3, data)

Any advice is appeciated.


cheers

Chris Akram



        [[alternative HTML version deleted]]

_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo

Re: intersections to individual polygons

Mon, 07/02/2018 - 22:44
Agustin Lobo writes:


>>> Given an SpatialPolygons DF with i.e. 2 intersecting polygons, A and B,
>>> is there a function that would split the original polygons into "onlyA", "onlyB" and
>>> "AintersectingB" polygons?


First convert the SpatialPolygons to spatstat 'owin' objects using as.owin.


If there are two polygons A and B, you can just use intersect.owin and setminus.owin to get what you want:

      AandB <- intersect.owin(A, B )

      AnotB <- setminus.owin(A, B )

      BnotA <- setminus.owin(B, A)


If there are several windows, make a list P containing the windows.

Then call

     Z <- kaleido(P)

     plot(Z, do.col=TRUE)

using the following code. The result Z is a tessellation. Extract the individual pieces as tiles(Z).


      kaleido <- function(P) {

           U <- union.owin(as.solist(P))

            V <- lapply(P, onesplit, U=U)

            Z <- Reduce(intersect.tess, V)

            return(Z)

      }

      onesplit <- function(X, U) tess(tiles=list(A=X, NotA=setminus.owin(U, X)), W=U)



Prof Adrian Baddeley DSc FAA

John Curtin Distinguished Professor

Department of Mathematics and Statistics

Curtin University, Perth, Western Australia



        [[alternative HTML version deleted]]

_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo

Re: intersections to individual polygons

Mon, 07/02/2018 - 20:07

Agustin:

I have attached a function that I think does what you want.
It returns either a tessellation (if you leave tess=TRUE) or
a list of owin objects otherwise.

I *think* that attachments with .R extensions make it through to the
list.  The attachment should at least get through to Agustin (who is
the person of central interest!).

It produces all non-empty intersections between sets of tiles
in the tessellation argument "ttt".  If singles=TRUE (the default)
these "intersections" include the tiles themselves.  Otherwise
at least two tiles are involved in each intersection.

E.g.

x1 <- allIntersections(te) # a tessellation
x2 <- allIntersections(te,tess=FALSE) # a list of windows
x3 <- allIntersections(te,singles=FALSE) # tiles themselves omitted.

Note that the code includes a work-around for an infelicity in
in intersect.owin().  This infelicity will very likely be fixed
in a future release of spatstat.

cheers,

Rolf

--
Technical Editor ANZJS
Department of Statistics
University of Auckland
Phone: +64-9-373-7599 ext. 88276

On 02/07/18 23:35, Agustin Lobo wrote:
> Mimicking your vignette:
>
> require(rgdal)
> require(raster)
> require(spatstat)
> require(rgeos)
> require(maptools)
>
> dzfile <-"https://www.dropbox.com/s/xxujcwqy3ec21sa/TDM1_DEM__04_S63W062_epsg32720.zip?dl=0"
> download.file(dzfile,"TDM1_DEM__04_S63W062_epsg32720.zip",method="wget")
> unzip("TDM1_DEM__04_S63W062_epsg32720.zip")
> v <- readOGR(dsn="TDM1_DEM__04_S63W062_epsg32720",
>               layer="TDM1_DEM__04_S63W062_epsg32720", stringsAsFactors = FALSE)
> plot(v)
> p <- slot(v, "polygons")
> p2 <- lapply(p, function(x) { SpatialPolygons(list(x)) })
> w <- lapply(p2, as.owin) #maptools required
>
> te <- tess(tiles=w)
> class(te)
> plot.tess(te,do.labels=TRUE)
>
> But this is still the original polygons,
> not the mosaic of polygon parts I'm looking for.
> Would I need to go doing all possible intersections and then
> adding the "A not B" and "B not A" parts for all possible pairs?
>
> Or is there a function in spatstat to convert "te" into a tesselation
> of adjacent polygons?
>
> Thanks
>
>
>
> On Thu, Jun 28, 2018 at 2:02 PM, Rolf Turner <[hidden email]> wrote:
>>
>> On 29/06/18 00:51, Agustin Lobo wrote:
>>
>>> Given an SpatialPolygons DF with i.e. 2 intersecting polygons, A and
>>> ,B,  is there a function
>>> that would split the original polygons into "onlyA", "onlyB" and
>>> "AintersectingB" polygons?
>>
>>
>> You can do this easily in spatstat by converting your polygons to "owin"
>> objects.
>>
>> cheers,
>>
>> Rolf Turner
>>
>> --
>> Technical Editor ANZJS
>> Department of Statistics
>> University of Auckland
>> Phone: +64-9-373-7599 ext. 88276
>>
>> _______________________________________________
>> R-sig-Geo mailing list
>> [hidden email]
>> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo

allIntersections.R (1K) Download Attachment

Re: Fwd: Spatiotemporal analysis r command

Mon, 07/02/2018 - 11:25
https://www.r-project.org/posting-guide.html

IHTH




Mauricio Zambrano-Bigiarini, PhD

=====================================
Department of Civil Engineering
Faculty of Engineering and Sciences
Universidad de La Frontera, Temuco, Chile
http://hzambran.github.io/
=====================================
mailto : [hidden email]
work-phone : +56 45 259 2812
=====================================
"Focus is about saying No" (Steve Jobs)
=====================================
Linux user #454569 -- Linux Mint user


On 1 July 2018 at 06:05, Yalemzewod Gelaw <[hidden email]> wrote:
> Hello everyone,
> I am planning to do a spatiotemporal distribution of hiv/tb coinfection and
> the effects of socio-demographic variables. I have a shapefile for 128
> districts.  For this ecological area  analysis I am planning to use
> Bayesian Poisson regression analysis.
> I here with looking for your help tutorial/notes/ r commands for the
> spatiotemporal analysis and Bayesian Poisson regression.
>
> Thank you for any help.
>
> *Regards, *
>
> Yalem
> --
> Best, Yalemzewod Assefa Gelaw (Yalem)
>
>         [[alternative HTML version deleted]]
>
> _______________________________________________
> R-sig-Geo mailing list
> [hidden email]
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
--
La información contenida en este correo electrónico y cualquier anexo o
respuesta relacionada, puede contener datos e información confidencial y no
puede ser usada o difundida por personas distintas a su(s) destinatario(s).
Si usted no es el destinatario de esta comunicación, le informamos que
cualquier divulgación, distribución o copia de esta información constituye
un delito conforme a la ley chilena. Si lo ha recibido por error, por favor
borre el mensaje y todos sus anexos y notifique al remitente.

_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo

Re: intersections to individual polygons

Mon, 07/02/2018 - 06:35
Mimicking your vignette:

require(rgdal)
require(raster)
require(spatstat)
require(rgeos)
require(maptools)

dzfile <-"https://www.dropbox.com/s/xxujcwqy3ec21sa/TDM1_DEM__04_S63W062_epsg32720.zip?dl=0"
download.file(dzfile,"TDM1_DEM__04_S63W062_epsg32720.zip",method="wget")
unzip("TDM1_DEM__04_S63W062_epsg32720.zip")
v <- readOGR(dsn="TDM1_DEM__04_S63W062_epsg32720",
             layer="TDM1_DEM__04_S63W062_epsg32720", stringsAsFactors = FALSE)
plot(v)
p <- slot(v, "polygons")
p2 <- lapply(p, function(x) { SpatialPolygons(list(x)) })
w <- lapply(p2, as.owin) #maptools required

te <- tess(tiles=w)
class(te)
plot.tess(te,do.labels=TRUE)

But this is still the original polygons,
not the mosaic of polygon parts I'm looking for.
Would I need to go doing all possible intersections and then
adding the "A not B" and "B not A" parts for all possible pairs?

Or is there a function in spatstat to convert "te" into a tesselation
of adjacent polygons?

Thanks



On Thu, Jun 28, 2018 at 2:02 PM, Rolf Turner <[hidden email]> wrote:
>
> On 29/06/18 00:51, Agustin Lobo wrote:
>
>> Given an SpatialPolygons DF with i.e. 2 intersecting polygons, A and
>> ,B,  is there a function
>> that would split the original polygons into "onlyA", "onlyB" and
>> "AintersectingB" polygons?
>
>
> You can do this easily in spatstat by converting your polygons to "owin"
> objects.
>
> cheers,
>
> Rolf Turner
>
> --
> Technical Editor ANZJS
> Department of Statistics
> University of Auckland
> Phone: +64-9-373-7599 ext. 88276
>
> _______________________________________________
> R-sig-Geo mailing list
> [hidden email]
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo

Fwd: Spatiotemporal analysis r command

Sun, 07/01/2018 - 05:05
Hello everyone,
I am planning to do a spatiotemporal distribution of hiv/tb coinfection and
the effects of socio-demographic variables. I have a shapefile for 128
districts.  For this ecological area  analysis I am planning to use
Bayesian Poisson regression analysis.
​I here with looking for your help tutorial/notes/ r commands for the
spatiotemporal analysis and Bayesian Poisson regression.​

​Thank you for any help.​

*Regards, *

Yalem
--
Best, Yalemzewod Assefa Gelaw (Yalem)

        [[alternative HTML version deleted]]

_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo

Re: Converting Map Grid Australia to Geographic coordinates

Sat, 06/30/2018 - 19:21
Brian,
You could hava a look at spTransform function in rgdal package.
Regards,
Zivan


> On 30 Jun 2018, at 15:55, Brian Williams <[hidden email]> wrote:
>
> Does anybody know of a package which can convert MGA coordinates to
> Geographic (lat/longs) ?
>
> Thanks
>
> Brian
>
>    [[alternative HTML version deleted]]
>
> _______________________________________________
> R-sig-Geo mailing list
> [hidden email]
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo

Converting Map Grid Australia to Geographic coordinates

Sat, 06/30/2018 - 14:55
Does anybody know of a package which can convert MGA coordinates to
Geographic (lat/longs) ?

Thanks

Brian

        [[alternative HTML version deleted]]

_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo

Re: intersections to individual polygons

Thu, 06/28/2018 - 08:02

On 29/06/18 00:51, Agustin Lobo wrote:

> Given an SpatialPolygons DF with i.e. 2 intersecting polygons, A and
> ,B,  is there a function
> that would split the original polygons into "onlyA", "onlyB" and
> "AintersectingB" polygons?

You can do this easily in spatstat by converting your polygons to "owin"
objects.

cheers,

Rolf Turner

--
Technical Editor ANZJS
Department of Statistics
University of Auckland
Phone: +64-9-373-7599 ext. 88276

_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo

intersections to individual polygons

Thu, 06/28/2018 - 07:51
Given an SpatialPolygons DF with i.e. 2 intersecting polygons, A and
,B,  is there a function
that would split the original polygons into "onlyA", "onlyB" and
"AintersectingB" polygons?

Thanks
Agus

_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo

Re: Fitting GWR for Panel Data in R

Wed, 06/27/2018 - 12:54
See

Fotheringham, A. S., Crespo, R., & Yao, J. (2015). Geographical and
Temporal Weighted Regression (GTWR). *Geographical Analysis*, 1–22.
https://doi.org/10.1111/gean.12071

Gollini, I., Lu, B., Charlton, M., Brunsdon, C., & Harris, P. (2015).
GWmodel: An R Package for Exploring Spatial Heterogeneity Using
Geographically Weighted Models. *Journal of Statistical Software*, *63*(17),
1–50. https://doi.org/10.1080/10095020.2014.917453
HTH, Dexter


On Wed, Jun 27, 2018 at 7:28 AM Arnold Salvacion via R-sig-Geo <
[hidden email]> wrote:

> Good day!
>
> Does anyone here have experience or code for fitting geographically
> weighted regression (GWR) for panel data. For example, I need to determine
> the effect of independent variable X1, X2, and X3 to Y which are all
> collected for T1 to Tn.  I have experience in fitting GWR for
> cross-sectional data but have no experience on doing such on panel data.
>
> Any suggestion is highly appreciated.
>
> Thanks in advance!
>
> Arnold Salvacion
>         [[alternative HTML version deleted]]
>
> _______________________________________________
> R-sig-Geo mailing list
> [hidden email]
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>
        [[alternative HTML version deleted]]

_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo

Re: Difference in area calculation between QGIS and R

Wed, 06/27/2018 - 11:31
 Thank you for you answer, Michael,
I am new to the use of R for spatial analyses.

I have always used QGIS for this kind of operations, but this time I need
to repeat the area calculation with thousands of polygons several times in
a loop, so I switched to R.

Actually, what my work needs is the most accurate measurement of such
areas. I am wondering if the best thing to do is to accept the results in R
and ignore those in QGIS.

Thanks,
Nick

2018-06-27 17:29 GMT+02:00 Michael Sumner <[hidden email]>:

> Raster uses a (discretized) cosine-of-latitude approximati (popular
> amongst longlat map makers).
>
> QGIS uses a project to local equal area projection method or maybe some
> other approach.
>
> There's lots and f options, all that matters is what your work needs.
>
> Cheers, Mike
>
> On Wed, 27 Jun 2018, 22:14 Suncus Etruscus, <[hidden email]>
> wrote:
>
>> Dear List,
>> I am trying to use R ("sp" and "raster" packages) to calculate the area of
>> several polygons (CRS of the shapefile EPSG: 4326  - WGS84).
>>
>> I used this line of code:
>>
>> shapefile_name$area_km2 <- area(shapefile_name)/1000000
>>
>> However, when I used QGIS to calculate the area of the same polygons
>> (through the "$area" function), I found there was a slight difference (in
>> every polygon).
>>
>> For example, for a polygon of about 30 000 km2, the area calculated in R
>> was 50 km2 smaller.
>>
>> What could be the cause?
>>
>> Thanks in advance,
>> N.
>>
>>         [[alternative HTML version deleted]]
>>
>> _______________________________________________
>> R-sig-Geo mailing list
>> [hidden email]
>> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>>
> --
> Dr. Michael Sumner
> Software and Database Engineer
> Australian Antarctic Division
> 203 Channel Highway
> <https://maps.google.com/?q=203+Channel+Highway+Kingston+Tasmania+7050+Australia&entry=gmail&source=g>
> Kingston Tasmania 7050 Australia
> <https://maps.google.com/?q=203+Channel+Highway+Kingston+Tasmania+7050+Australia&entry=gmail&source=g>
>
>
        [[alternative HTML version deleted]]

_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo

Re: Difference in area calculation between QGIS and R

Wed, 06/27/2018 - 11:07

On 06/27/2018 05:29 PM, Michael Sumner wrote:
> Raster uses a (discretized) cosine-of-latitude approximati (popular amongst
> longlat map makers).

As far as I can tell, raster branches into the libgeographic code
(copied into the package sources, although Karney is not mentioned as
one of the package copyright holders), using close to exact geodesic /
ellipsoidal computation.

>
> QGIS uses a project to local equal area projection method or maybe some
> other approach.
>
> There's lots and f options, all that matters is what your work needs.
>
> Cheers, Mike
>
> On Wed, 27 Jun 2018, 22:14 Suncus Etruscus, <[hidden email]>
> wrote:
>
>> Dear List,
>> I am trying to use R ("sp" and "raster" packages) to calculate the area of
>> several polygons (CRS of the shapefile EPSG: 4326  - WGS84).
>>
>> I used this line of code:
>>
>> shapefile_name$area_km2 <- area(shapefile_name)/1000000
>>
>> However, when I used QGIS to calculate the area of the same polygons
>> (through the "$area" function), I found there was a slight difference (in
>> every polygon).
>>
>> For example, for a polygon of about 30 000 km2, the area calculated in R
>> was 50 km2 smaller.
>>
>> What could be the cause?
>>
>> Thanks in advance,
>> N.
>>
>>         [[alternative HTML version deleted]]
>>
>> _______________________________________________
>> R-sig-Geo mailing list
>> [hidden email]
>> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>>
--
Edzer Pebesma
Institute for Geoinformatics
Heisenbergstrasse 2, 48151 Muenster, Germany
Phone: +49 251 8333081

_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo

Pages