[R-sig-Geo] correct use of autoKrige and gstat with Digital Elevation Models

Stefano Sofia @tef@no@@ofi@ @ending from regione@m@rche@it
Fri Aug 17 15:56:10 CEST 2018


Hi Zivan.
You are right. Yes, I did it:

Inverse_Distance_rainfall <- gstat(formula=Cumulata~ExtractedElevationValues, data=df_prec, locations=newgrid)

My problem is when I try to put the output in a raster. How can I do it?
Inverse_Distance_rainfall is a list, but within it I am not able to find the name to use within raster().
Any suggestion?

Thank you
Stefano




         (oo)
--oOO--( )--OOo----------------
Stefano Sofia PhD
Area Meteorologica e  Area nivologica - Centro Funzionale
Servizio Protezione Civile - Regione Marche
Via del Colle Ameno 5
60126 Torrette di Ancona, Ancona
Uff: 071 806 7743
E-mail: stefano.sofia using regione.marche.it
---Oo---------oO----------------
________________________________________
Da: Zivan Karaman [zivan.karaman using gmail.com]
Inviato: venerdì 17 agosto 2018 14.50
A: Stefano Sofia
Oggetto: Re: [R-sig-Geo] correct use of autoKrige and gstat with Digital Elevation Models

Have you tried:
Inverse_Distance_rainfall <- gstat(formula=Cumulata~ExtractedElevationValues, data=df_prec, locations=newgrid)
formula, data and locations are NOT the default first 3 arguments to gstat function - see help(gstat).
HTH
Zivan

> Le 16 août 2018 à 09:21, Stefano Sofia <stefano.sofia using regione.marche.it> a écrit :
>
> Dear R-sig-geo list users,
> I am dealing with rainfall isohyet maps, I am using a DEM in order to get more accurate values.
> In particular I am trying to compare Kriging method (through autoKrige) with the Inverse Distance method (through gstat).
> Kriging code works fine, I think that I am not able to make a correct use of gstat.
>
> In both cases,
> - I load my rainfall cumulates in a data frame called df_prec:
>
> Codice_sensore Cumulata Long_Cent Lat_Cent Quota
> 3134 177 13.22888 42.99361  1352
> 3135 98 13.19046 42.85843  1496
> 3136 40 12.81444 43.23840  1005
> 3137 201 13.30681 42.84369  1036
> 3138 145 13.32638 42.90215   980
> 3139 45 13.11180 42.79578  1575
> 3140 69 13.18538 42.91439  1800
> 3234 129 13.44506 42.74417   960
>
> - I load the shape file and the DEM:
>
> ## LOAD SHAPEFILE
> Shape_area5 <- readOGR(dsn="area5.shp", layer="area5")
> Shape_area5_projection <- "+init=epsg:32633 +proj=utm +zone=33 +datum=WGS84 +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0"
> Shape_area5 <- spTransform(Shape_area5, CRS(Shape_area5_projection))
> ## LOADING DEM
> ita_DEM <- getData('alt', country='ITA', mask=FALSE)
> ita_DEM <- projectRaster(ita_DEM, crs = CRS("+init=epsg:32633"))
>
> - I put the coordinates in the correct way, I evaluate the extracted elevation values and I create a new grid with these new values:
>
> ##
> coordinates(df_prec) <- ~Long_Cent+Lat_Cent
> proj4string(df_prec) <- CRS("+init=epsg:4326")
> df_prec <- spTransform(df_prec, CRS("+init=epsg:32633"))
> ## EXTRACT THE ELEVATION VALUES TO MY POINTS; extract EXPECTS A RASTER INPUT
> df_prec$ExtractedElevationValues <- extract(x=ita_DEM, y=df_prec)
> ## CREATE A NEW GRID
> newgrid <- as(ita_DEM, "SpatialGridDataFrame")
> names(newgrid) <- "ExtractedElevationValues"
>
> - when I apply Kriging, everything works fine
> ## KRIGING
> RegressionKriging_rainfall <- autoKrige(Cumulata~ExtractedElevationValues, df_prec, newgrid)
> prec_rast <- raster(RegressionKriging_rainfall$krige_output)
>
> - if I use gstat in the standard way
>
> Inverse_Distance_rainfall <- gstat(formula=Cumulata~1, data=df_prec, locations=newgrid)
>
> everything goes fine, but this code obviously does not take into consideration the elevation points. When I try to use the inverse distance method like Kriging:
> Inverse_Distance_rainfall <- gstat(formula=Cumulata~ExtractedElevationValues, df_prec, newgrid)
>
> I get the following error:
>
> Error in (function (classes, fdef, mtable)  :
>  unable to find an inherited method for function ‘proj4string’ for signature ‘"NULL"’
>  Calls: my_inversedistanceweighted_DEM -> gstat -> proj4string -> <Anonymous>
>
>
> Where is my mistake?
> Could please somebody help me and give me a hint for a correst use of gstat?
> Thank you for your attention and your help
> Stefano
>
>
>         (oo)
> --oOO--( )--OOo----------------
> Stefano Sofia PhD
> Area Meteorologica e  Area nivologica - Centro Funzionale
> Servizio Protezione Civile - Regione Marche
> Via del Colle Ameno 5
> 60126 Torrette di Ancona, Ancona
> Uff: 071 806 7743
> E-mail: stefano.sofia using regione.marche.it
> ---Oo---------oO----------------
>
> ________________________________
>
> AVVISO IMPORTANTE: Questo messaggio di posta elettronica può contenere informazioni confidenziali, pertanto è destinato solo a persone autorizzate alla ricezione. I messaggi di posta elettronica per i client di Regione Marche possono contenere informazioni confidenziali e con privilegi legali. Se non si è il destinatario specificato, non leggere, copiare, inoltrare o archiviare questo messaggio. Se si è ricevuto questo messaggio per errore, inoltrarlo al mittente ed eliminarlo completamente dal sistema del proprio computer. Ai sensi dell’art. 6 della DGR n. 1394/2008 si segnala che, in caso di necessità ed urgenza, la risposta al presente messaggio di posta elettronica può essere visionata da persone estranee al destinatario.
> IMPORTANT NOTICE: This e-mail message is intended to be received only by persons entitled to receive the confidential information it may contain. E-mail messages to clients of Regione Marche may contain information that is confidential and legally privileged. Please do not read, copy, forward, or store this message unless you are an intended recipient of it. If you have received this message in error, please forward it to the sender and delete it completely from your computer system.
>
> --
> Questo messaggio  stato analizzato da Libra ESVA ed  risultato non infetto.
> This message was scanned by Libra ESVA and is believed to be clean.
>
>
>    [[alternative HTML version deleted]]
>
> _______________________________________________
> R-sig-Geo mailing list
> R-sig-Geo using r-project.org
>  https://urlsand.esvalabs.com/?u=https%3A%2F%2Fstat.ethz.ch%2Fmailman%2Flistinfo%2Fr-sig-geo&e=52342f8a&h=f7d4a649&f=y&p=y

--

Questo messaggio  stato analizzato con Libra ESVA ed  risultato non infetto.


________________________________

AVVISO IMPORTANTE: Questo messaggio di posta elettronica può contenere informazioni confidenziali, pertanto è destinato solo a persone autorizzate alla ricezione. I messaggi di posta elettronica per i client di Regione Marche possono contenere informazioni confidenziali e con privilegi legali. Se non si è il destinatario specificato, non leggere, copiare, inoltrare o archiviare questo messaggio. Se si è ricevuto questo messaggio per errore, inoltrarlo al mittente ed eliminarlo completamente dal sistema del proprio computer. Ai sensi dell’art. 6 della DGR n. 1394/2008 si segnala che, in caso di necessità ed urgenza, la risposta al presente messaggio di posta elettronica può essere visionata da persone estranee al destinatario.
IMPORTANT NOTICE: This e-mail message is intended to be received only by persons entitled to receive the confidential information it may contain. E-mail messages to clients of Regione Marche may contain information that is confidential and legally privileged. Please do not read, copy, forward, or store this message unless you are an intended recipient of it. If you have received this message in error, please forward it to the sender and delete it completely from your computer system.

--
Questo messaggio  stato analizzato da Libra ESVA ed  risultato non infetto.
This message was scanned by Libra ESVA and is believed to be clean.



More information about the R-sig-Geo mailing list