knitr::opts_chunk$set(echo = TRUE, cache = FALSE)
Look at Figure \@ref(fig:testing), it shows the output of a magic function.
library(TestPackage) a <- add_2(seq(1:10)) plot(a)
Adding 2 to 4 we get r add_2(4)
.
knitr::kable(head(cars), caption = "A table of the cars package" )
Look at the table \@ref(tab:table).
Thanks to @bookdown, @rmarkdown, @knitr
Xie knows what he (?) is doing [-@knitr].
add_3(5)
library(tidyverse) x <- rnorm(1e6) - 5e3 knitr::kable(head(x))
knitr::kable(head(x + 2))
selsecond(mtcars)
library(binford) data(LRB) knitr::kable(head(LRB))
harran <- read.table("../data/Sites_HarranPlain.csv", sep = ",", header=TRUE) str(harran) library(sp) coordinates(harran) <- ~X+Y proj4string(harran) <- CRS("+init=epsg:4326") str(harran)
library(raster) ##srtm <- getData("SRTM", lon=38, lat=37) srtm <- raster("srtm_44_05.tif") plot(srtm) points(harran) srtm <- crop(srtm, extent(harran)+1) plot(srtm) srtm <- projectRaster(srtm, crs = CRS("+init=epsg:32637")) srtm2 <- aggregate(srtm, fact = 2) writeRaster(srtm2, "data/dem.tif", overwrite = TRUE)
harran <- spTransform(harran, CRSobj = CRS("+init=epsg:32637")) library(spatstat) harran_ppp <- ppp(x = harran@coords[,1], y = harran@coords[,2], window = owin(xrange = harran@bbox[1,], yrange = c(min(harran@coords[,2]), min(harran@coords[,2]+52000)))) plot(harran)
anyDuplicated(harran_ppp) harran <- unique(harran_ppp) plot(harran_ppp) harran_ppp_nn <- nndist(harran_ppp) str(harran_ppp_nn) hist(harran_ppp_nn)
harran_kde <- density.ppp(x = harran_ppp, sigma = mean(harran_ppp_nn)) plot(harran_kde)
library(raster) dem <- raster("../data/dem.tif") im_dem <- as.im(as.image.SpatialGridDataFrame(as(dem, "SpatialGridDataFrame"))) harran_rhohat <- rhohat(object = harran_ppp, covariate = im_dem, bw = 200 ) plot(harran_rhohat) str(harran_rhohat) rho_dem <- predict(harran_rhohat) plot(rho_dem) diff_rho <- harran_kde - rho_dem
create random points with rpoispp function that have the same intensity like our point pattern.
set.seed(123) harran_poispp2 <- rpoispp(lambda = harran_ppp$n/area.owin(harran_ppp$window), win = harran_ppp$window) set.seed(123) harran_poispp3 <- rpoispp(intensity(harran_ppp), win=Window(harran_ppp)) set.seed(123) harran_poispp4 <- rpoispp(ex = harran_ppp) plot(harran_ppp) points(harran_poispp2, col = "red") points(harran_poispp3, col = "blue") points(harran_poispp4, col = "green")
harran_g <- Gest(harran_ppp) str(harran_g) plot(harran_g) harran_ge <- envelope(harran_ppp, fun = "Gest", nsim = 999) plot(harran_ge) harran_f <- Fest(harran_ppp) harran_fe <- envelope(harran_ppp, fun = "Fest", nsim = 99) harran_k <- Kest(harran_ppp) harran_ke <- envelope(harran_ppp, fun = "Kest", nsim = 99) plot(harran_fe) plot(harran_ke)
harran_gi <- Ginhom(harran_ppp, lambda = predict(harran_rhohat)) #par(mfrow = c(1,2)) #plot(harran_gi, xlim = c(0,6000)) #plot(harran_g, xlim = c(0,6000))
load("../data/Precipitation.RData") test <- data.frame(test) library(sp) coordinates(test) <- ~lon+lat proj4string(test) <- CRS("+init=epsg:4326") test2 <- spTransform(test, CRS("+init=epsg:32634")) plot(test2) library(raster) ##srtm <- getData("SRTM", lon=mean(coordinates(test)[,1]), lat=mean(coordinates(test)[,2])) srtm <- raster("../data/srtm_41_05.tif") srtm <- crop(srtm, extent(test)+1) plot(srtm) srtm3 <- projectRaster(srtm, crs=CRS("+init=epsg:32634")) srtm3 <- aggregate(srtm3, fact = 3) library(gstat)
rain_idw <- idw(mean_r~1, test2, as(srtm3, "SpatialGridDataFrame")) rain_idw_cv <- krige.cv(mean_r~1, test2) head(rain_idw_cv) plot(rain_idw)
plot(variogram(mean_r~1, loc = test2, cloud = TRUE)) va <- variogram(mean_r~1, loc = test2, cutoff = 80000, width = 5000) plot(va) plot(va, vgm(8e+04, "Exp", 5e+04, 5000)) fva <- fit.variogram(va, vgm(8e+04, "Exp", 5e+04, 5000), fit.method = 7) ## ordinary least squares; default plot(va, fva) rain_krige <- krige(mean_r~1, test2, as(srtm3, "SpatialGridDataFrame"), fva) rain_krige_cv <- krige.cv(mean_r~1, test2, fva) head(rain_krige_cv) plot(rain_krige) rain_krige2 <- brick(rain_krige) plot(rain_krige2)
library(automap) rainautkrige <- automap::autoKrige(mean_r~1, test2, as(srtm3, "SpatialGridDataFrame")) str(rainautkrige) plot(rainautkrige) tmp <- rainautkrige[[1]] tmp <- raster(tmp) tmp2 <- raster(rain_krige) tmp3 <- tmp2-tmp plot(tmp3) ## or using Reduce function (that needs a list) plot(Reduce("-", list(tmp2,tmp)))
names(srtm3) <- "altitude" test3 <- test2[!is.na(test2$altitude),] rainkwed <- automap::autoKrige(mean_r~altitude, test3, as(srtm3, "SpatialGridDataFrame")) plot(rainkwed)
see Wolfgang's script.
linmod <- lm(mean_r~altitude, test3) linmod plot(variogram(linmod$residuals~1, loc = test3, cloud = TRUE))
see Wolfgang's script.
test3$x <- coordinates(test3)[,1]
library(lattice); library(gstat) v <- variogram(mean_r ~ 1, test3) xyplot( x = gamma/1e8 ~ dist, data = v, pch = 3, type = 'b', lwd = 2, panel = function(x, y, ...) { for (i in 1:500) { test3$random = sample(test3$mean_r) v = variogram(random ~ 1, test3) llines(x = v$dist, y = v$gamma/1e8, col = rgb(.5,.5,.5,.2) ) } panel.xyplot(x, y, ...) }, xlab = 'distance', ylab = 'semivariance/1e8' )
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.