Mengubah peta raster tipe habitat secara acak?

12

Saya memiliki raster tipe habitat untuk wilayah tertentu di Skotlandia. Saya perlu membuat skenario habitat masa depan dengan perubahan habitat untuk menilai kelayakan populasi spesies burung.

Misalnya, di masa depan mungkin ada 10% lebih banyak kehutanan di daerah tersebut. Saya ingin mengubah peta saat ini dengan menambahkan secara acak kehutanan dalam blok dengan ukuran tertentu. Sejauh ini saya berpikir untuk memilih titik acak dari raster yang mengidentifikasi area di mana kehutanan dapat terjadi dan menumbuhkan blok berukuran benar menggunakan semacam automata seluler.

Apakah ini sepertinya cara terbaik untuk melakukan ini? Apakah ada metode yang lebih baik?

Jika ini cara terbaik yang tersedia, bagaimana saya bisa melakukan ini, lebih disukai, R? (Saat ini saya sedang melihat fungsi rpoint di "spatstat" bersama dengan paket CellularAutomata)

Saya juga memiliki akses ke GRASS, QGis dan ArcMap 10 jika ada cara yang lebih sederhana di antara mereka.

Matt Geary
sumber
Apakah Anda sudah melihat rasterpaketnya? Itu punya banyak alat untuk bekerja dengan data raster (noo, rly?).
Roman Luštrik
Terima kasih, Roman. Ya, ini seharusnya memberi saya alat untuk membaca dan memanipulasi peta dasar saya.
Matt Geary

Jawaban:

18

Pernahkah Anda berpikir untuk menggunakan rantai Markov ? Ini secara efektif merupakan "otomat seluler seluler probabilistik," sehingga memasok keacakan yang diinginkan. Alih-alih meresepkan generasi baru dalam hal tetangga lokal dari generasi yang ada, itu menentukan distribusi probabilitas untuk generasi baru. Distribusi itu dapat diperkirakan dari, katakanlah, urutan waktu gambar dari area yang sama atau serupa.

Secara intuitif, model ini mengatakan bahwa sel tidak harus membuat transisi dari berhutan ke tidak berhutan (atau sebaliknya ), tetapi kemungkinan bahwa hal itu akan membuat transisi tergantung pada tutupan lahan segera di sekitarnya. Ia dapat menangani beberapa kelas tutupan, konfigurasi kompleks lingkungan, dan bahkan digeneralisasikan untuk "mengingat" sejarah evolusi tutupan lahan baru-baru ini.

Transisi dapat diimplementasikan menggunakan pernyataan Peta Aljabar, yang membuat metode ini praktis dalam GIS berbasis raster, bahkan yang tanpa akses langsung atau cepat ke data tingkat sel. Menggunakan R membuatnya lebih mudah.

Sebagai contoh, pertimbangkan konfigurasi awal ini hanya dengan dua kelas, putih dan hitam:

Grid tutupan lahan

Untuk menggambarkan apa yang bisa terjadi, saya membuat model parameter (tidak didasarkan pada data apa pun) di mana transisi ke hitam terjadi dengan probabilitas 1 - q ^ k di mana k adalah jumlah rata-rata sel hitam dalam lingkungan 3 oleh 3 (k = 0, 1/9, 2/9, ..., 1). Ketika q kecil atau sebagian besar lingkungan sudah hitam, sel baru akan hitam. Berikut adalah empat simulasi independen dari generasi kesepuluh untuk lima nilai q mulai dari 0,25 hingga 0,05:

Daftar hasil

Jelas model ini memiliki banyak karakteristik CA tetapi juga mencakup efek acak yang berguna untuk mengeksplorasi hasil alternatif.


Kode

Berikut ini mengimplementasikan simulasi di R.

#
# Make a transition from state `x` using a kernel having `k.ft` as
# its Fourier transform.
#
transition <- function(x, k.ft, q=0.1) {
  k <- zapsmall(Re(fft(k.ft * fft(x), inverse=TRUE))) / length(x)
  matrix(runif(length(k)) > q^k, nrow=nrow(k))
}
#
# Create the zeroth generation and the fft of a transition kernel.
#
n.row <- 2^7 # FFT is best with powers of 2
n.col <- 2^7
kernel <- matrix(0, nrow=n.row, ncol=n.col)
kernel[1:3, 1:3] <- 1/9
kernel.f <- fft(kernel)

set.seed(17)
x <- matrix(sample(c(0,1), n.row*n.col, replace=TRUE, prob=c(599, 1)), n.row)
#
# Prepare to run multiple simulations.
#
y.list <- list()
parameters <- c(.25, .2, .15, .1, .05)
#
# Perform and benchmark the simulations.
#
i <- 0
system.time({
  for (q in parameters) {
    y <- x
    for (generation in 1:10) {
      y <- transition(y, kernel.f, q)
    }
    y.list[[i <- i+1]] <- y
  }
})
#
# Display the results.
#    
par(mfrow=c(1,length(parameters)))
invisible(sapply(1:length(parameters), 
       function(i) image(y.list[[i]], 
                         col=c("White", "Black"),
                         main=parameters[i])))
whuber
sumber
+1 Sangat menarik. Jika Anda memiliki data tutupan lahan bersejarah untuk area tertentu, apakah mungkin untuk mendapatkan q dan / atau k?
Kirk Kuykendall
2
@Kirk Ya, tapi saya tidak akan merekomendasikannya: model yang saya gunakan untuk ilustrasi terlalu sederhana. Tetapi Anda dapat memperoleh sesuatu yang lebih baik: dengan melihat frekuensi transisi empiris dari setiap konfigurasi lingkungan yang telah terjadi, Anda dapat membuat model evolusi masa depan yang transisinya secara statistik meniru evolusi masa lalu. Jika frekuensi transisi secara homogen spasial dan jika masa depan terus bertindak seperti masa lalu, menjalankan beberapa simulasi ini dapat memberikan gambaran yang jelas tentang apa yang mungkin terjadi di masa depan.
whuber
Terima kasih, ini sepertinya melakukan apa yang saya butuhkan. Apakah mungkin untuk menetapkan batas proporsi area yang berubah?
Matt Geary
@ Mat Ya, setidaknya dalam arti probabilistik. Teori ini menjelaskan bagaimana rantai Markov mencapai campuran yang stabil secara asimptot dari masing-masing negara bagian. Ini adalah keseimbangan dinamis: pada setiap generasi banyak sel mungkin berubah, tetapi hasil akhirnya adalah menjaga proporsinya dalam grid tetap sama (hingga penyimpangan peluang kecil).
whuber
1
Saya seorang programmer R yang mengerikan. Saya dapat membagikan kode Mathematica yang saya gunakan; dengan fungsi yang diterapkan R, itu harus port dengan baik. Anda memerlukan kernel, aturan transisi, dan prosedur untuk menerapkannya ke array 2D 0/1. Jadi: kernel = ConstantArray[1/3^2, {3,3}]untuk kernel; transitionRule [k_] := With[{q = 0.1}, Boole[RandomReal[{0, 1}] > q^k]]untuk aturan; dan next[a_, kernel_, f_] := Map[f, ListConvolve[kernel, a, {1, 1}, 0], {2}]untuk menerapkannya ke array a . Misalnya, untuk merencanakan empat generasi dari awal , gunakan ArrayPlot /@ NestList[next[#, kernel, transitionRule] &, start, 3].
whuber