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: 36 min 34 sec ago

Re: spatstat - kppm function

Tue, 06/06/2017 - 18:43
On 07/06/17 09:57, Fernanda De Bastiani wrote:

>Dear Sir/Madam,
>I am trying to use function kppm() of spatstat package with anisotropic
>covariance option as the following:
>#-------------------------------------------
>    kppm(redwood, ~x, "LGCP", statistic="pcf",
>covmodel=list(model="matern",
>nu=0.3, Aniso = matrix(nc=2, c(1.5, 3, -3, 4))))
>#-------------------------------------------
>but it gives me an error message:
>Error in rfInit(model = p, x = x, y = y, z = z, T = T, grid = grid,
>distances = distances,  :
>    'RFcov' :  'Aniso' : not of the required size: (2, 2) instead of
>(undetermined, 1)
>Please, could anyone help me to understand in which conditions nc=2 will
>work?
The kppm() function is designed for models whose pair correlation function
is isotropic.
It cannot be used to fit an anisotropic model.

(Extending this to anisotropic models is conceivable but it would require
new theory as well as completely new software).

I¹m sorry about the opaque error message. I¹ll try to find a way to detect
anisotropic models and give a more meaningful message from the spatstat
code.

Adrian Baddeley

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

spatstat - kppm function

Tue, 06/06/2017 - 16:57
Dear Sir/Madam,

I am trying to use function kppm() of spatstat package with anisotropic
covariance option as the following:
#-------------------------------------------
  kppm(redwood, ~x, "LGCP", statistic="pcf",  covmodel=list(model="matern",
nu=0.3, Aniso = matrix(nc=2, c(1.5, 3, -3, 4))))
#-------------------------------------------
but it gives me an error message:

Error in rfInit(model = p, x = x, y = y, z = z, T = T, grid = grid,
distances = distances,  :
  'RFcov' :  'Aniso' : not of the required size: (2, 2) instead of
(undetermined, 1)

Please, could anyone help me to understand in which conditions nc=2 will
work?

Many thanks,
Fernanda

        [[alternative HTML version deleted]]

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

Post Question

Sat, 06/03/2017 - 08:44
I would like to post this question:
I am trying to figure out how to reassign small polygons to its nearest neighbor in R (similar to gdal_sieve.py function).
I created a spatial polygon file from a classified gridded dataset using the grid2polygons function in R. I would now like to take that polygon file and merge any polygons containing less than 6 pixels with its nearest largest neighbor. In the image link (https://i.stack.imgur.com/6BLvP.jpg), you can see an example of small polygons circled in blue. I would like those merged with the larger surround polygon. Alternatively, if there is a way to create a minimum threshold of 6 clumped pixels prior to the conversion to polygon stage, any suggestions for that method would be great as well.
There is a function in Gdal: gdal_sieve.py that can do exactly this (http://www.gdal.org/gdal_sieve.html) if you would like to get an idea of exactly what I would like to do. Does anyone know of a way to do this in R?
        [[alternative HTML version deleted]]

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

Re: Problem in ENFA: infinite or missing values

Thu, 06/01/2017 - 22:53
Hey Perry,

It is not clear to me where the problem lies without more
details/explanations. As Mike said, can you come up with a reproducible
example? If you can't make it minimal, see if you can send me the data
privately, I'll try to investigate on my side.

Cheers,
Mathieu.


On 05/19/2017 02:51 AM, Perry Beasley-hall wrote:
> Hello everyone,
> I'm trying to do an ENFA using the package adehabitatHS. I'm following the example code and am creating the 'pr' file using:
>
> pr <- slot(count.points(species, map), "data")[,1]
>
> Where 'species' is a .csv file with two columns of latitude and longitude, and 'map' is a SpatialGridDataFrame file of 23 .asc inputs. For reference, I'm importing my .asc files using readGDAL and then using the cbind command to put them into the 'map' file.
>
> I encounter the following warning message when trying to create 'pr':
>
> 1: In points2grid(points, tolerance, round) :
>   grid has empty column/rows in dimension 1
> 2: In points2grid(points, tolerance, round) :
>   grid has empty column/rows in dimension 2
>
> I assumed this would be alright, but I'm getting the subsequent error when I attempt the enfa command:
>
> enfa1 <- enfa(pc, pr, scannf = FALSE)
>
> Error in eigen(Se) : infinite or missing values in 'x'
>
> Some searching online has told me that this might mean not all of my rows in the 'species' .csv correspond to variables in the 'map' file, but I don't think this is the case for my dataset.
>
> Can anyone help me with this issue? Many thanks in advance.
> ____
>
> Perry Beasley-Hall
> PhD Candidate
> School of Life and Environmental Sciences
>
> THE UNIVERSITY OF SYDNEY
> Room 329, Edgeworth David Building A11
> Sydney NSW 2006, Australia
> E-mail: [hidden email]
> Phone: +61 420 494 584
>
> [[alternative HTML version deleted]]
>
> _______________________________________________
> R-sig-Geo mailing list
> [hidden email]
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>
--

Mathieu Basille

[hidden email] | http://ase-research.org/basille
+1 954-577-6314 | University of Florida FLREC

  « Le tout est de tout dire, et je manque de mots
  Et je manque de temps, et je manque d'audace. »
  — Paul Éluard

This message is signed to guarantee its authenticity.
For a true private correspondence, use my public key
to encrypt your messages:

  http://mathieu.basille.net/pub.asc

Learn more: http://mzl.la/1BsOGiZ

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

Re: Google Earth Engine?

Tue, 05/30/2017 - 12:45

On 5/26/17, 4:02 PM, "[hidden email] on behalf of Barry
Rowlingson" <[hidden email] on behalf of
[hidden email]> wrote:

>On Fri, May 26, 2017 at 11:34 PM, Andy Bunn <[hidden email]> wrote:
>> Does anybody out there interface with the google earth engine from R?
>>I'm too old a dog to learn python. -Andy
>
>Too old? Never! See: https://www.xkcd.com/353/


Perfect! I had forgotten that one. Well, I think I might have found a new
project. The GEE seems like the perfect tool for Landsat products. Thanks,
A



>
>Given that the other supported option is Javascript....
>
>I suspect a solution using one of the R-Python interfaces
>("reticulate" perhaps) might be the best solution, but you still might
>have to learn python to construct the analysis jobs.
>
>Barry
>
>
>
>>         [[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

Re: How to fit a pure spatial variogram on a spatio-temporal empirical one

Tue, 05/30/2017 - 07:26
Thank you very much Ben, the code works perfectly.
I have 2 years time series, one for each station, of daily mean average pm10 which I have detrended wrt mean and standard deviation by using land use coefficients.
The problem I am trying to face is to distinguish between the spatial and the temporal components of the correlation among time series, the former due to distance, the latter mostly due to seasonality.
Does the pooled variogram obtained with the code you posted here isolates the temporal component of the correlation and only assume constant covariance structure over time or is it influenced by seasonality?

> Il giorno 30 mag 2017, alle ore 10:48, Dr. Benedikt Gräler <[hidden email]> ha scritto:
>
> Dear Carlo,
>
> the code below is a bit of a hack, but does what you are asking for. The classes "gstatVariogram" and "StVariogram" have slight different design and so do the functions fit.variogram and fit.StVariogram. Note that spVv is now a pooled variogram across all time steps of your dataset treating each time slice as an independent copy of the same pure spatial process (i.e. strong temporal autocorrelation might influence your estimation).
>
> HTH,
>
> Ben
>
>
> library(gstat)
> data("vv")
> plot(vv)
>
> spaceOnly <- vv$timelag == 0
>
> spVv <- cbind(vv[spaceOnly,],
>              data.frame(dir.hor=rep(0, sum(spaceOnly)),
>                         dir.ver=rep(0, sum(spaceOnly))))
>
> # drop empty (NA) first row
> spVv <- spVv[-1, ]
>
> # manually re-class
> class(spVv) <- c("gstatVariogram","data.frame")
>
> plot(spVv)
>
> fitSpVgm <- fit.variogram(spVv, vgm(30, "Exp", 150, 10))
> plot(spVv, fitSpVgm)
>
>
>
>> On 29/05/2017 20:13, Carlo Cavalieri wrote:
>> Hi, I am using the R package GSTAT to make a spatio-temporal interpolation for my thesis and I wanted to know if it was possible to obtain the pure spatial empirical variogram from the spatio-temporal so that I can use it to fit a pure spatial variogram, for example exponential.
>> Unfortunately fit.variogram only accepts objects output of variogram, not of variogramST. One possible solution could be to extract tlag=0 from the StVariogram and convert the output to class variogramModel, but I have no idea on how to do this.
>> I look for a way to do this because fit the spatial variogram for each day separately is not a good idea given the small number of observation stations.
>> One way is definitely possible since the authors of the paper "Spatio-Temporal Interpolation using gstat” managed to compare the results of pure spatial and spatio-temporal interpolation (that is what I want to do): Below a quotation from that paper.
>> "For comparison with classical approaches, we interpolate across Germany iteratively for each single day using all available data for variogram estimation. The purely spatial empirical variogram can directly be obtained from the empirical spatio-temporal variogram, by fixing the temporal lag at 0 separation. From the same set of variogram models as investigated for the spatio-temporal models, the exponential model (partial sill: 66.5, range: 224 km, nugget: 13.5) is the best suited based on the optimisation criterion.”
>> Does anyone have any idea?
>> Thank you
>> Carlo
>>    [[alternative HTML version deleted]]
>> _______________________________________________
>> R-sig-Geo mailing list
>> [hidden email]
>> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>
> --
> Dr. Benedikt Gräler
> 52°North Initiative for Geospatial Open Source Software GmbH
> Martin-Luther-King-Weg 24
> 48155 Muenster, Germany
>
> E-Mail: [hidden email]
> Fon: +49-(0)-251/396371-39
> Fax: +49-(0)-251/396371-11
>
> http://52north.org/
> Twitter: @FiveTwoN
>
> General Managers: Dr. Albert Remke, Dr. Andreas Wytzisk
> Local Court Muenster HRB 10849
> <b_graeler.vcf>
> _______________________________________________
> 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

Re: Deal with multiple factorlevel in one grid square

Tue, 05/30/2017 - 07:07
maybe fasterize? Only on GitHub, and requires sf.

(Mixing raster and sf is piece-meal, but very doable).

https://github.com/ecohealthalliance/fasterize

On Tue, 30 May 2017, 21:55 Miriam Püts <[hidden email]> wrote:

> Hi all,
>
> @Pat: yes, that is almost what I need. If I would do it in ArcGIS I would
> use Polygon to raster (with a second raster as reference) and then choose
> majority.  there a way to do it in R? If I use raster::extract() I would
> extract data from the raster object for the locations of other spatial
> data, but I need spatial data for the raster object...
>
>
> Cheers,
> Miriam
>
> <°)))>< >°)))>< >°)))>< >°)))><
>
> Miriam Püts
> Marine Lebende Ressourcen/ Marine Living Resources
> Thünen-Institut für Seefischerei/ Thünen Institute of Sea fisheries
> Palmaille 9
> 22767 Hamburg (Germany)
>
> Tel:  +49 40 38905-105
> Mail: [hidden email]
>
> ----- Ursprüngliche Mail -----
> Von: "Alexander Herr" <[hidden email]>
> An: "patrick schratz" <[hidden email]>, [hidden email],
> "miriam puets" <[hidden email]>
> CC: [hidden email]
> Gesendet: Donnerstag, 25. Mai 2017 02:28:40
> Betreff: RE: [R-sig-Geo] Deal with multiple factorlevel in one grid square
>
> Hiya,
>
> Do you want a raster with attribute information (ie several attributes
> from your  shape file)? You can achieve this with a workflow that uses
> raster::extract() to assigne unique IDs of polygon to raster cells
> create a RAT for you raster
> assign to the unique IDs of raster RAT the attributes of your polygon
> unique IDs (many to one relationship) with something like merge() or
> functions in libraries like diplyr or data.table
>
> Cheers
> Herry
>
> -----Original Message-----
> From: R-sig-Geo [mailto:[hidden email]] On Behalf Of
> Patrick Schratz
> Sent: Wednesday, 24 May 2017 10:48 PM
> To: Michael Sumner <[hidden email]>; Miriam Püts <
> [hidden email]>
> Cc: [hidden email]
> Subject: Re: [R-sig-Geo] Deal with multiple factorlevel in one grid square
>
> Do you want to do what is called “zonal statistics” in ArcGIS with
> “majority" option?
>
> You may check out this SO question and try the “mode” function within
> raster::extract() - maybe it does what you need.
>
> (I’m also unsure if I understand the question correctly)
>
> Cheers, Pat
>
> PhD Student at Department of Geography - GIScience group
> Friedrich-Schiller-University Jena, Germany
> Tel.: +49-3641-9-48973
> Web: https://pat-s.github.io
>
> On 24. May 2017, 12:35 +0200, Miriam Püts <[hidden email]>,
> wrote:
> > Hi Mike,
> >
> > I will try to explain it a bit more in detail. Maybe it is easier to
> understand if I start from the end. In the end I would like to have an
> ASCII file to read into Ecospace, which has the same extensions and
> coordinates as other files I already created. This ASCII should contain
> information about the sediment type within each predefined cell. To create
> this ASCII file I have a shape file with the polygons representing the
> sediment type and my grid which I applied to other variables to have the
> same extend. Now I would like to create a grid containing the information
> on sediment. Here, per grid cell the sediment type which covers the most of
> the cell should be defined and connected with the coordinates for the grid
> cell.
> >
> > I hope this makes it more clear...
> >
> >
> > <°)))>< >°)))>< >°)))>< >°)))><
> >
> > Miriam Püts
> > Marine Lebende Ressourcen/ Marine Living Resources Thünen-Institut für
> > Seefischerei/ Thünen Institute of Sea fisheries Palmaille 9
> > 22767 Hamburg (Germany)
> >
> > Tel: +49 40 38905-105
> > Mail: [hidden email]
> >
> >
> > Von: "Michael Sumner" <[hidden email]
> > An: "Miriam Püts" <[hidden email]>, [hidden email]
> > Gesendet: Mittwoch, 24. Mai 2017 11:55:14
> > Betreff: Re: [R-sig-Geo] Deal with multiple factorlevel in one grid
> > square
> >
> >
> >
> > raster::extract(grid, poly, weights = TRUE) is a start. It returns a
> list which is painful to deal with at first, but can be collected into one
> data frame for standard summarizing.
> >
> > I'm still a bit confused about whether you want an estimate of a cell
> overlap in a polygon or something else.
> >
> > Cheers, Mike
> > On Wed, 24 May 2017, 17:12 Miriam Püts < [ mailto:
> [hidden email] | [hidden email] ] > wrote:
> >
> >
> > Hi everyone,
> >
> > I have the following problem: I have a personalized grid and a shape
> file with polygons representing sediment types. Now I would like to apply
> this grid to the Polygons to identify the sediment type most common within
> each grid. I tried it with rasterize, but here I can only chose last or
> first. Die you have any suggestions how I might get the sediment type for
> each raster cell?
> >
> > Thank you for your help!
> >
> > <°)))>< >°)))>< >°)))>< >°)))><
> >
> > Miriam Püts
> > Marine Lebende Ressourcen/ Marine Living Resources Thünen-Institut für
> > Seefischerei/ Thünen Institute of Sea fisheries Palmaille 9
> > 22767 Hamburg (Germany)
> >
> > Tel: +49 40 38905-105
> > Mail: [ mailto:[hidden email] | [hidden email] ]
> >
> > _______________________________________________
> > R-sig-Geo mailing list
> > [ mailto:[hidden email] | [hidden email] ] [
> > https://stat.ethz.ch/mailman/listinfo/r-sig-geo |
> > https://stat.ethz.ch/mailman/listinfo/r-sig-geo ]
> >
> >
> >
> > --
> > Dr. Michael Sumner
> > Software and Database Engineer
> > Australian Antarctic Division
> > 203 Channel Highway
> > Kingston Tasmania 7050 Australia
> >
> >
> > [[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
> --
Dr. Michael Sumner
Software and Database Engineer
Australian Antarctic Division
203 Channel Highway
Kingston Tasmania 7050 Australia

        [[alternative HTML version deleted]]

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

Re: Deal with multiple factorlevel in one grid square

Tue, 05/30/2017 - 06:55
Hi all,

@Pat: yes, that is almost what I need. If I would do it in ArcGIS I would use Polygon to raster (with a second raster as reference) and then choose majority.  there a way to do it in R? If I use raster::extract() I would extract data from the raster object for the locations of other spatial data, but I need spatial data for the raster object...


Cheers,
Miriam

<°)))>< >°)))>< >°)))>< >°)))><

Miriam Püts
Marine Lebende Ressourcen/ Marine Living Resources
Thünen-Institut für Seefischerei/ Thünen Institute of Sea fisheries
Palmaille 9
22767 Hamburg (Germany)

Tel:  +49 40 38905-105
Mail: [hidden email]

----- Ursprüngliche Mail -----
Von: "Alexander Herr" <[hidden email]>
An: "patrick schratz" <[hidden email]>, [hidden email], "miriam puets" <[hidden email]>
CC: [hidden email]
Gesendet: Donnerstag, 25. Mai 2017 02:28:40
Betreff: RE: [R-sig-Geo] Deal with multiple factorlevel in one grid square

Hiya,

Do you want a raster with attribute information (ie several attributes from your  shape file)? You can achieve this with a workflow that uses
raster::extract() to assigne unique IDs of polygon to raster cells
create a RAT for you raster
assign to the unique IDs of raster RAT the attributes of your polygon unique IDs (many to one relationship) with something like merge() or functions in libraries like diplyr or data.table

Cheers
Herry

-----Original Message-----
From: R-sig-Geo [mailto:[hidden email]] On Behalf Of Patrick Schratz
Sent: Wednesday, 24 May 2017 10:48 PM
To: Michael Sumner <[hidden email]>; Miriam Püts <[hidden email]>
Cc: [hidden email]
Subject: Re: [R-sig-Geo] Deal with multiple factorlevel in one grid square

Do you want to do what is called “zonal statistics” in ArcGIS with “majority" option?

You may check out this SO question and try the “mode” function within raster::extract() - maybe it does what you need.

(I’m also unsure if I understand the question correctly)

Cheers, Pat

PhD Student at Department of Geography - GIScience group Friedrich-Schiller-University Jena, Germany
Tel.: +49-3641-9-48973
Web: https://pat-s.github.io

On 24. May 2017, 12:35 +0200, Miriam Püts <[hidden email]>, wrote:
> Hi Mike,
>
> I will try to explain it a bit more in detail. Maybe it is easier to understand if I start from the end. In the end I would like to have an ASCII file to read into Ecospace, which has the same extensions and coordinates as other files I already created. This ASCII should contain information about the sediment type within each predefined cell. To create this ASCII file I have a shape file with the polygons representing the sediment type and my grid which I applied to other variables to have the same extend. Now I would like to create a grid containing the information on sediment. Here, per grid cell the sediment type which covers the most of the cell should be defined and connected with the coordinates for the grid cell.
>
> I hope this makes it more clear...
>
>
> <°)))>< >°)))>< >°)))>< >°)))><
>
> Miriam Püts
> Marine Lebende Ressourcen/ Marine Living Resources Thünen-Institut für
> Seefischerei/ Thünen Institute of Sea fisheries Palmaille 9
> 22767 Hamburg (Germany)
>
> Tel: +49 40 38905-105
> Mail: [hidden email]
>
>
> Von: "Michael Sumner" <[hidden email]
> An: "Miriam Püts" <[hidden email]>, [hidden email]
> Gesendet: Mittwoch, 24. Mai 2017 11:55:14
> Betreff: Re: [R-sig-Geo] Deal with multiple factorlevel in one grid
> square
>
>
>
> raster::extract(grid, poly, weights = TRUE) is a start. It returns a list which is painful to deal with at first, but can be collected into one data frame for standard summarizing.
>
> I'm still a bit confused about whether you want an estimate of a cell overlap in a polygon or something else.
>
> Cheers, Mike
> On Wed, 24 May 2017, 17:12 Miriam Püts < [ mailto:[hidden email] | [hidden email] ] > wrote:
>
>
> Hi everyone,
>
> I have the following problem: I have a personalized grid and a shape file with polygons representing sediment types. Now I would like to apply this grid to the Polygons to identify the sediment type most common within each grid. I tried it with rasterize, but here I can only chose last or first. Die you have any suggestions how I might get the sediment type for each raster cell?
>
> Thank you for your help!
>
> <°)))>< >°)))>< >°)))>< >°)))><
>
> Miriam Püts
> Marine Lebende Ressourcen/ Marine Living Resources Thünen-Institut für
> Seefischerei/ Thünen Institute of Sea fisheries Palmaille 9
> 22767 Hamburg (Germany)
>
> Tel: +49 40 38905-105
> Mail: [ mailto:[hidden email] | [hidden email] ]
>
> _______________________________________________
> R-sig-Geo mailing list
> [ mailto:[hidden email] | [hidden email] ] [
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo |
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo ]
>
>
>
> --
> Dr. Michael Sumner
> Software and Database Engineer
> Australian Antarctic Division
> 203 Channel Highway
> Kingston Tasmania 7050 Australia
>
>
> [[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

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

Re: How to fit a pure spatial variogram on a spatio-temporal empirical one

Tue, 05/30/2017 - 03:48
Dear Carlo,

the code below is a bit of a hack, but does what you are asking for. The
classes "gstatVariogram" and "StVariogram" have slight different design
and so do the functions fit.variogram and fit.StVariogram. Note that
spVv is now a pooled variogram across all time steps of your dataset
treating each time slice as an independent copy of the same pure spatial
process (i.e. strong temporal autocorrelation might influence your
estimation).

HTH,

  Ben


library(gstat)
data("vv")
plot(vv)

spaceOnly <- vv$timelag == 0

spVv <- cbind(vv[spaceOnly,],
               data.frame(dir.hor=rep(0, sum(spaceOnly)),
                          dir.ver=rep(0, sum(spaceOnly))))

# drop empty (NA) first row
spVv <- spVv[-1, ]

# manually re-class
class(spVv) <- c("gstatVariogram","data.frame")

plot(spVv)

fitSpVgm <- fit.variogram(spVv, vgm(30, "Exp", 150, 10))
plot(spVv, fitSpVgm)



On 29/05/2017 20:13, Carlo Cavalieri wrote:
> Hi, I am using the R package GSTAT to make a spatio-temporal interpolation for my thesis and I wanted to know if it was possible to obtain the pure spatial empirical variogram from the spatio-temporal so that I can use it to fit a pure spatial variogram, for example exponential.
> Unfortunately fit.variogram only accepts objects output of variogram, not of variogramST. One possible solution could be to extract tlag=0 from the StVariogram and convert the output to class variogramModel, but I have no idea on how to do this.
> I look for a way to do this because fit the spatial variogram for each day separately is not a good idea given the small number of observation stations.
>
> One way is definitely possible since the authors of the paper "Spatio-Temporal Interpolation using gstat” managed to compare the results of pure spatial and spatio-temporal interpolation (that is what I want to do): Below a quotation from that paper.
>
> "For comparison with classical approaches, we interpolate across Germany iteratively for each single day using all available data for variogram estimation. The purely spatial empirical variogram can directly be obtained from the empirical spatio-temporal variogram, by fixing the temporal lag at 0 separation. From the same set of variogram models as investigated for the spatio-temporal models, the exponential model (partial sill: 66.5, range: 224 km, nugget: 13.5) is the best suited based on the optimisation criterion.”
>
> Does anyone have any idea?
> Thank you
> Carlo
> [[alternative HTML version deleted]]
>
> _______________________________________________
> R-sig-Geo mailing list
> [hidden email]
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
> --
Dr. Benedikt Gräler
52°North Initiative for Geospatial Open Source Software GmbH
Martin-Luther-King-Weg 24
48155 Muenster, Germany

E-Mail: [hidden email]
Fon: +49-(0)-251/396371-39
Fax: +49-(0)-251/396371-11

http://52north.org/
Twitter: @FiveTwoN

General Managers: Dr. Albert Remke, Dr. Andreas Wytzisk
Local Court Muenster HRB 10849

_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo
b_graeler.vcf (435 bytes) Download Attachment

How to fit a pure spatial variogram on a spatio-temporal empirical one

Mon, 05/29/2017 - 13:13
Hi, I am using the R package GSTAT to make a spatio-temporal interpolation for my thesis and I wanted to know if it was possible to obtain the pure spatial empirical variogram from the spatio-temporal so that I can use it to fit a pure spatial variogram, for example exponential.
Unfortunately fit.variogram only accepts objects output of variogram, not of variogramST. One possible solution could be to extract tlag=0 from the StVariogram and convert the output to class variogramModel, but I have no idea on how to do this.
I look for a way to do this because fit the spatial variogram for each day separately is not a good idea given the small number of observation stations.

One way is definitely possible since the authors of the paper "Spatio-Temporal Interpolation using gstat” managed to compare the results of pure spatial and spatio-temporal interpolation (that is what I want to do): Below a quotation from that paper.

"For comparison with classical approaches, we interpolate across Germany iteratively for each single day using all available data for variogram estimation. The purely spatial empirical variogram can directly be obtained from the empirical spatio-temporal variogram, by fixing the temporal lag at 0 separation. From the same set of variogram models as investigated for the spatio-temporal models, the exponential model (partial sill: 66.5, range: 224 km, nugget: 13.5) is the best suited based on the optimisation criterion.”

Does anyone have any idea?
Thank you
Carlo
        [[alternative HTML version deleted]]

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

Re: Google Earth Engine?

Fri, 05/26/2017 - 18:02
On Fri, May 26, 2017 at 11:34 PM, Andy Bunn <[hidden email]> wrote:
> Does anybody out there interface with the google earth engine from R? I'm too old a dog to learn python. -Andy

Too old? Never! See: https://www.xkcd.com/353/

Given that the other supported option is Javascript....

I suspect a solution using one of the R-Python interfaces
("reticulate" perhaps) might be the best solution, but you still might
have to learn python to construct the analysis jobs.

Barry



>         [[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

Google Earth Engine?

Fri, 05/26/2017 - 17:34
Does anybody out there interface with the google earth engine from R? I'm too old a dog to learn python. -Andy

        [[alternative HTML version deleted]]

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

Re: release of RQGIS 1.0.0

Fri, 05/26/2017 - 17:25
Wonderful work! Thank you for your contributions!

On 5/26/17, 1:06 AM, "R-sig-Geo on behalf of "Jannes Münchow""
<[hidden email] on behalf of [hidden email]> wrote:

>Today we are proud to announce a major RQGIS release (already on CRAN).
>We have completely rewritten RQGIS now using the reticulate package to
>access the QGIS Python API. Overall, RQGIS is now even more
>user-friendly. Find more information on new features, what has changed
>and how to use RQGIS on following pages:
>- https://jannesm.wordpress.com/category/r/
>- https://github.com/jannes-m/RQGIS/releases
>- https://github.com/jannes-m/RQGIS
>
>_______________________________________________
>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

release of RQGIS 1.0.0

Fri, 05/26/2017 - 03:06
Today we are proud to announce a major RQGIS release (already on CRAN). We have completely rewritten RQGIS now using the reticulate package to access the QGIS Python API. Overall, RQGIS is now even more user-friendly. Find more information on new features, what has changed and how to use RQGIS on following pages:
- https://jannesm.wordpress.com/category/r/
- https://github.com/jannes-m/RQGIS/releases
- https://github.com/jannes-m/RQGIS

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

Re: Easiest way to get a a subset from MODIS

Wed, 05/24/2017 - 23:34
Answering my own question.

I was able to do what I wanted using the MODIS package. The runGdal function does what I want. Thanks to those who coded it!

-A

From: Andy Bunn <[hidden email]<mailto:[hidden email]>>
Date: Wednesday, May 24, 2017 at 12:22 PM
To: R-sig-Geo <[hidden email]<mailto:[hidden email]>>
Subject: Easiest way to get a a subset from MODIS

Hi all, It's been years since I've used any MODIS data. Looking around I see there are many R packages that purport to interface nicely with MODIS servers and deliver you a subset of data.

I have an EarthData login and can (and have) download the tile I want via http://search.earthdata.nasa.gov but getting the hdf, converting to tif, clipping etc. is tedious and I'm wondering what the best way to do this in R is. So,  let's imagine I would like to get a class(raster) object of the "250m_16_days_NDVI" band from "MYD13Q1" product for a 10x10km area around the point c(68.75,161.42) for the period of nearest July 1, 2015. What package and approach is the easiest in R?

Is it something like this?

library(MODISTools)
dat2get <- data.frame(lat=68.75,long=161.42,
                      start.date=as.POSIXlt("2015-07-01"),
                      end.date=as.POSIXlt("2015-07-01"),
                      ID=1)

MODISSubsets(LoadDat = dat2get,
             Products = "MYD13Q1",
             Bands = "250m_16_days_NDVI",
             Size = c(10,10),
             StartDate = TRUE)

Regards, Andy



        [[alternative HTML version deleted]]

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

Re: Deal with multiple factorlevel in one grid square

Wed, 05/24/2017 - 19:28
Hiya,

Do you want a raster with attribute information (ie several attributes from your  shape file)? You can achieve this with a workflow that uses
raster::extract() to assigne unique IDs of polygon to raster cells
create a RAT for you raster
assign to the unique IDs of raster RAT the attributes of your polygon unique IDs (many to one relationship) with something like merge() or functions in libraries like diplyr or data.table

Cheers
Herry

-----Original Message-----
From: R-sig-Geo [mailto:[hidden email]] On Behalf Of Patrick Schratz
Sent: Wednesday, 24 May 2017 10:48 PM
To: Michael Sumner <[hidden email]>; Miriam Püts <[hidden email]>
Cc: [hidden email]
Subject: Re: [R-sig-Geo] Deal with multiple factorlevel in one grid square

Do you want to do what is called “zonal statistics” in ArcGIS with “majority" option?

You may check out this SO question and try the “mode” function within raster::extract() - maybe it does what you need.

(I’m also unsure if I understand the question correctly)

Cheers, Pat

PhD Student at Department of Geography - GIScience group Friedrich-Schiller-University Jena, Germany
Tel.: +49-3641-9-48973
Web: https://pat-s.github.io

On 24. May 2017, 12:35 +0200, Miriam Püts <[hidden email]>, wrote:
> Hi Mike,
>
> I will try to explain it a bit more in detail. Maybe it is easier to understand if I start from the end. In the end I would like to have an ASCII file to read into Ecospace, which has the same extensions and coordinates as other files I already created. This ASCII should contain information about the sediment type within each predefined cell. To create this ASCII file I have a shape file with the polygons representing the sediment type and my grid which I applied to other variables to have the same extend. Now I would like to create a grid containing the information on sediment. Here, per grid cell the sediment type which covers the most of the cell should be defined and connected with the coordinates for the grid cell.
>
> I hope this makes it more clear...
>
>
> <°)))>< >°)))>< >°)))>< >°)))><
>
> Miriam Püts
> Marine Lebende Ressourcen/ Marine Living Resources Thünen-Institut für
> Seefischerei/ Thünen Institute of Sea fisheries Palmaille 9
> 22767 Hamburg (Germany)
>
> Tel: +49 40 38905-105
> Mail: [hidden email]
>
>
> Von: "Michael Sumner" <[hidden email]
> An: "Miriam Püts" <[hidden email]>, [hidden email]
> Gesendet: Mittwoch, 24. Mai 2017 11:55:14
> Betreff: Re: [R-sig-Geo] Deal with multiple factorlevel in one grid
> square
>
>
>
> raster::extract(grid, poly, weights = TRUE) is a start. It returns a list which is painful to deal with at first, but can be collected into one data frame for standard summarizing.
>
> I'm still a bit confused about whether you want an estimate of a cell overlap in a polygon or something else.
>
> Cheers, Mike
> On Wed, 24 May 2017, 17:12 Miriam Püts < [ mailto:[hidden email] | [hidden email] ] > wrote:
>
>
> Hi everyone,
>
> I have the following problem: I have a personalized grid and a shape file with polygons representing sediment types. Now I would like to apply this grid to the Polygons to identify the sediment type most common within each grid. I tried it with rasterize, but here I can only chose last or first. Die you have any suggestions how I might get the sediment type for each raster cell?
>
> Thank you for your help!
>
> <°)))>< >°)))>< >°)))>< >°)))><
>
> Miriam Püts
> Marine Lebende Ressourcen/ Marine Living Resources Thünen-Institut für
> Seefischerei/ Thünen Institute of Sea fisheries Palmaille 9
> 22767 Hamburg (Germany)
>
> Tel: +49 40 38905-105
> Mail: [ mailto:[hidden email] | [hidden email] ]
>
> _______________________________________________
> R-sig-Geo mailing list
> [ mailto:[hidden email] | [hidden email] ] [
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo |
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo ]
>
>
>
> --
> Dr. Michael Sumner
> Software and Database Engineer
> Australian Antarctic Division
> 203 Channel Highway
> Kingston Tasmania 7050 Australia
>
>
> [[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
_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo

Re: Recs for R tools that subset/tile raster bricks?

Wed, 05/24/2017 - 18:18
Thanks for your response, Mike! The files are ENVI images (*.dat). The
output of print() on one is copied below. The other images have slightly
different row/col sizes but the same number of layers.

Thanks for the tip on the vignettes. It looks like "Writing functions with
the ”raster” package" is the one you're referring to? blockSize() looks
helpful as I'm not sure what size would be appropriate to crop to.

--

class       : RasterBrick
dimensions  : 24390, 2032, 49560480, 373  (nrow, ncol, ncell, nlayers)
resolution  : 1, 1  (x, y)
extent      : 520554.6, 522586.6, 4766639, 4791029  (xmin, xmax, ymin, ymax)
coord. ref. : +proj=utm +zone=11 +datum=NAD83 +units=m +no_defs
+ellps=GRS80 +towgs84=0,0,0
data source : E:\RCEW_CoReg_2014\183421\183421_Georegistered
names       :
Warp..Band.8.ang20140919t183421_corr_v1c_img...381.360260.Nanometers.,
Warp..Band.9.ang20140919t183421_corr_v1c_img...386.368927.Nanometers.,
Warp..Band.10.ang20140919t183421_corr_v1c_img...391.377594.Nanometers.,
Warp..Band.11.ang20140919t183421_corr_v1c_img...396.386292.Nanometers.,
Warp..Band.12.ang20140919t183421_corr_v1c_img...401.394958.Nanometers.,
Warp..Band.13.ang20140919t183421_corr_v1c_img...406.403625.Nanometers.,
Warp..Band.14.ang20140919t183421_corr_v1c_img...411.412292.Nanometers.,
Warp..Band.15.ang20140919t183421_corr_v1c_img...416.420959.Nanometers.,
Warp..Band.16.ang20140919t183421_corr_v1c_img...421.429626.Nanometers.,
Warp..Band.17.ang20140919t183421_corr_v1c_img...426.438293.Nanometers.,
Warp..Band.18.ang20140919t183421_corr_v1c_img...431.446960.Nanometers.,
Warp..Band.19.ang20140919t183421_corr_v1c_img...436.455627.Nanometers.,
Warp..Band.20.ang20140919t183421_corr_v1c_img...441.464294.Nanometers.,
Warp..Band.21.ang20140919t183421_corr_v1c_img...446.... <truncated>

Thanks,
Megan

On Wed, May 24, 2017 at 6:58 PM, Michael Sumner <[hidden email]> wrote:

> Do also see the raster vignettes, there is one on performance and using
> line- or block- based chunks for large jobs.
>
> Cheers, Mike
>
>
> On Thu, 25 May 2017, 08:32 Michael Sumner <[hidden email]> wrote:
>
>> Try crop() (or more abstractly extract(raster, cells) once you've
>> determined which cells x-y are needed).
>>
>> Getting the most efficient way to do this depends on the source of the
>> brick. It's not possible to know without knowing details on your setup. The
>> file format, method of creation and dimensions of the brick and sizes of
>> your subset query can all affect the best ways to go.
>>
>> Happy to help but it might take a bit of conversation to figure out
>> what's best. If you can share the output of print on the brick anf the
>> extent of a subset that is a good start.
>>
>> Cheers, Mike
>>
>> On Thu, 25 May 2017, 05:50 Megan Maloney <[hidden email]> wrote:
>>
>>> Hello all,
>>>
>>> Does anyone know of a good R tool to subset/tile raster bricks?
>>>
>>> I'd like to break a raster brick into smaller pieces to process. I'm
>>> finding functions to subset or tile a single raster layer, but not for
>>> multi-layer raster bricks.
>>>
>>> The brick input is hyperspectral flightlines - very large, 1 m^2
>>> resampling
>>> covering a watershed, 340 layers. My code currently takes in a raster
>>> brick, extracts values per cell, applies beta coefficients, and outputs a
>>> nitrogen value to a raster.
>>>
>>> I don't need the subsets/tiles to be a particular size or number, just
>>> smaller. I was warned against dropping a whole flightline into R. I'm
>>> new to
>>> R, so I'm not sure what limitations I'm running into except processing
>>> time.
>>> All advice is welcome.
>>>
>>> Thank you very much for your help.
>>> Megan
>>>
>>>         [[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
>> Kingston Tasmania 7050 Australia
>>
>> --
> Dr. Michael Sumner
> Software and Database Engineer
> Australian Antarctic Division
> 203 Channel Highway
> Kingston Tasmania 7050 Australia
>
>
        [[alternative HTML version deleted]]

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

Re: Recs for R tools that subset/tile raster bricks?

Wed, 05/24/2017 - 17:58
Do also see the raster vignettes, there is one on performance and using
line- or block- based chunks for large jobs.

Cheers, Mike

On Thu, 25 May 2017, 08:32 Michael Sumner <[hidden email]> wrote:

> Try crop() (or more abstractly extract(raster, cells) once you've
> determined which cells x-y are needed).
>
> Getting the most efficient way to do this depends on the source of the
> brick. It's not possible to know without knowing details on your setup. The
> file format, method of creation and dimensions of the brick and sizes of
> your subset query can all affect the best ways to go.
>
> Happy to help but it might take a bit of conversation to figure out what's
> best. If you can share the output of print on the brick anf the extent of a
> subset that is a good start.
>
> Cheers, Mike
>
> On Thu, 25 May 2017, 05:50 Megan Maloney <[hidden email]> wrote:
>
>> Hello all,
>>
>> Does anyone know of a good R tool to subset/tile raster bricks?
>>
>> I'd like to break a raster brick into smaller pieces to process. I'm
>> finding functions to subset or tile a single raster layer, but not for
>> multi-layer raster bricks.
>>
>> The brick input is hyperspectral flightlines - very large, 1 m^2
>> resampling
>> covering a watershed, 340 layers. My code currently takes in a raster
>> brick, extracts values per cell, applies beta coefficients, and outputs a
>> nitrogen value to a raster.
>>
>> I don't need the subsets/tiles to be a particular size or number, just
>> smaller. I was warned against dropping a whole flightline into R. I'm new
>> to
>> R, so I'm not sure what limitations I'm running into except processing
>> time.
>> All advice is welcome.
>>
>> Thank you very much for your help.
>> Megan
>>
>>         [[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
> Kingston Tasmania 7050 Australia
>
> -- Dr. Michael Sumner
Software and Database Engineer
Australian Antarctic Division
203 Channel Highway
Kingston Tasmania 7050 Australia

        [[alternative HTML version deleted]]

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

Re: Recs for R tools that subset/tile raster bricks?

Wed, 05/24/2017 - 17:33
Try crop() (or more abstractly extract(raster, cells) once you've
determined which cells x-y are needed).

Getting the most efficient way to do this depends on the source of the
brick. It's not possible to know without knowing details on your setup. The
file format, method of creation and dimensions of the brick and sizes of
your subset query can all affect the best ways to go.

Happy to help but it might take a bit of conversation to figure out what's
best. If you can share the output of print on the brick anf the extent of a
subset that is a good start.

Cheers, Mike

On Thu, 25 May 2017, 05:50 Megan Maloney <[hidden email]> wrote:

> Hello all,
>
> Does anyone know of a good R tool to subset/tile raster bricks?
>
> I'd like to break a raster brick into smaller pieces to process. I'm
> finding functions to subset or tile a single raster layer, but not for
> multi-layer raster bricks.
>
> The brick input is hyperspectral flightlines - very large, 1 m^2 resampling
> covering a watershed, 340 layers. My code currently takes in a raster
> brick, extracts values per cell, applies beta coefficients, and outputs a
> nitrogen value to a raster.
>
> I don't need the subsets/tiles to be a particular size or number, just
> smaller. I was warned against dropping a whole flightline into R. I'm new
> to
> R, so I'm not sure what limitations I'm running into except processing
> time.
> All advice is welcome.
>
> Thank you very much for your help.
> Megan
>
>         [[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
Kingston Tasmania 7050 Australia

        [[alternative HTML version deleted]]

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

Recs for R tools that subset/tile raster bricks?

Wed, 05/24/2017 - 14:49
Hello all,

Does anyone know of a good R tool to subset/tile raster bricks?

I'd like to break a raster brick into smaller pieces to process. I'm
finding functions to subset or tile a single raster layer, but not for
multi-layer raster bricks.

The brick input is hyperspectral flightlines - very large, 1 m^2 resampling
covering a watershed, 340 layers. My code currently takes in a raster
brick, extracts values per cell, applies beta coefficients, and outputs a
nitrogen value to a raster.

I don't need the subsets/tiles to be a particular size or number, just
smaller. I was warned against dropping a whole flightline into R. I'm new to
R, so I'm not sure what limitations I'm running into except processing time.
All advice is welcome.

Thank you very much for your help.
Megan

        [[alternative HTML version deleted]]

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

Pages