correct use of autoKrige and gstat with Digital Elevation Models

classic Classic list List threaded Threaded
2 messages Options
Reply | Threaded
Open this post in threaded view
|

correct use of autoKrige and gstat with Digital Elevation Models

Stefano Sofia
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: [hidden email]
---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
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo
Reply | Threaded
Open this post in threaded view
|

Re: correct use of autoKrige and gstat with Digital Elevation Models

Stefano Sofia
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: [hidden email]
---Oo---------oO----------------
________________________________________
Da: Zivan Karaman [[hidden email]]
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 <[hidden email]> 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: [hidden email]
> ---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
> [hidden email]
>  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.

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