Contoh kriging biasa selangkah demi selangkah?

8

Saya telah mengikuti tutorial online untuk sprigal kriging dengan keduanya geoRdan gstat(dan juga automap). Saya dapat melakukan kriging spasial dan saya memahami konsep utama di baliknya. Saya tahu cara membuat semivariogram, cara menyesuaikan model dengan itu dan bagaimana melakukan kriging biasa.

Yang tidak saya mengerti adalah bagaimana bobot dari nilai yang diukur di sekitarnya ditentukan. Saya tahu mereka berasal dari semivariogram dan tergantung pada jarak dari lokasi prediksi dan pada pengaturan spasial dari titik-titik yang diukur. Tapi bagaimana caranya?

Adakah yang bisa membuat model kriging (non-bayesian) biasa dengan 3 titik pengukuran acak dan 1 lokasi prediksi? Itu akan mencerahkan.

Pigna
sumber
1
hanya karena penasaran, mengapa Anda tidak ingin melihat jawaban Bayesian? Itu membuat banyak hal lebih sederhana ketika Anda berurusan dengan Proses Gaussian.
DeltaIV
@DeltaIV karena pertama saya ingin belajar dengan cara yang sering. Statistik Bayesian masih mendung bagi saya
Pigna
1
" Yang tidak saya mengerti adalah bagaimana bobot dari nilai yang diukur di sekitarnya ditentukan. " Jika ada yang tertarik, saya mengirim jawaban di GIS SE dengan sebuah contoh tentang cara menghitungnya ( gis.stackexchange.com/questions/270274/… ). Tapi jawabannya di sini sudah bagus!
Andre Silva

Jawaban:

9

Pertama saya akan menjelaskan kriging biasa dengan tiga poin secara matematis. Asumsikan kita memiliki bidang acak intrinsik stasioner.

Kriging Biasa

Kami mencoba memprediksi nilainya Z(x0) menggunakan nilai yang diketahui Z=(Z(x1),Z(x2),Z(x3)) Prediksi yang kita inginkan adalah dari bentuk

Z^(x0)=λTZ
dimana λ=(λ1,λ2,λ3)adalah bobot interpolasi. Kami mengasumsikan nilai rata-rata konstanμ. Untuk mendapatkan hasil yang tidak bias, kami perbaikiλ1+λ2+λ3=1. Kami kemudian mendapatkan masalah berikut:
minE(Z(X0)λTZ)2s.t.λT1=1.
Menggunakan metode pengali Lagrange, kami memperoleh persamaan:
j=13λjγ(xsaya-xj)+m=γ(xsaya-x0),saya=1,2,3,
j=13λj=1,
dimana m adalah pengali lagrange dan γadalah variogram (semi). Dari ini, kita dapat mengamati beberapa hal:
  • Bobot tidak tergantung pada nilai rata-rata μ.
  • Bobot tidak tergantung pada nilai Zsama sekali. Hanya pada koordinat (dalam kasus isotropik pada jarak saja)
  • Setiap bobot tergantung pada lokasi semua titik lainnya.

Perilaku bobot yang tepat sulit dilihat hanya dari persamaannya, tetapi orang dapat dengan sangat kasar mengatakan:

  • Lebih lanjut intinya adalah dari x0, semakin rendah beratnya ("lebih lanjut" sehubungan dengan poin lain).
  • Namun, dekat dengan titik lain juga menurunkan berat badan.
  • Hasilnya sangat tergantung pada bentuk, jangkauan, dan, khususnya, efek nugget dari variogram. Akan sangat mencerahkan untuk mempertimbangkan krigingR dengan hanya dua poin dan lihat bagaimana hasilnya berubah dengan pengaturan variogram yang berbeda.

Namun saya akan fokus pada lokasi titik di pesawat. Saya menulis fungsi R kecil ini yang mengambil poin dari[0,1]2 dan plot bobot kriging (untuk fungsi kovarian eksponensial dengan nol nugget).

library(geoR)

# Plots prediction weights for kriging in the window [0,1]x[0,1] with the prediction point (0.5,0.5)
drawWeights <- function(x,y){
 df <- data.frame(x=x,y=y, values = rep(1,length(x)))
  data <- as.geodata(df, coords.col = 1:2, data.col = 3)

  wls <- variofit(bin1,ini=c(1,0.5),fix.nugget=T)
  weights <- round(as.numeric(krweights(data$coords,c(0.5,0.5),krige.control(obj.mod=wls, type="ok"))),3)

  plot(data$coords, xlim=c(0,1),  ylim=c(0,1))
  segments(rep(0.5,length(x)), rep(0.5,length(x)),x, y, lty=3 )
  text((x+0.5)/2,(y+0.5)/2,labels=weights)
}

Anda dapat bermain dengannya menggunakan clickpppfungsi spatstat :

library(spatstat)
points <- clickppp()
drawWeights(points$x,points$y)

Berikut ini beberapa contohnya

Poin berjarak sama dari x0 dan dari satu sama lain

deg <- seq(0,2*pi,length.out=4)
deg <- head(deg,length(deg)-1)
x <- 0.5*as.numeric(lapply(deg, cos)) + 0.5
y <- 0.5*as.numeric(lapply(deg, sin)) + 0.5
drawWeights(x,y)

masukkan deskripsi gambar di sini

Poin yang dekat satu sama lain akan berbagi bobot

deg <- c(0,0.1,pi)
x <- 0.5*as.numeric(lapply(deg, cos)) + 0.5
y <- 0.5*as.numeric(lapply(deg, sin)) + 0.5
drawWeights(x,y)

masukkan deskripsi gambar di sini

Titik terdekat "mencuri" bobot

deg <- seq(0,2*pi,length.out=4)
deg <- head(deg,length(deg)-1)
x <- c(0.6,0.5*as.numeric(lapply(deg, cos)) + 0.5)
y <- c(0.6,0.5*as.numeric(lapply(deg, sin)) + 0.5)
drawWeights(x,y)

https://i.imgur.com/MeuPvFT.png

Dimungkinkan untuk mendapatkan bobot negatif

masukkan deskripsi gambar di sini

Semoga ini memberi Anda merasakan bagaimana bobot bekerja.

Dahn
sumber