Poisson process
Kiranna Indy 5052241014
2026-05-26
3.2.
3.1. Exponential Interarrival
Simulasi lintasan (trajectory) proses Poisson dengan parameter laju λ=2 dan visualisasi
#specifying parameters
lambda<- 2
njumps<- 20#defining states
N<- 0:njumps#setting time as vector
time<- c()#setting initial value for time
time[1]<- 0#specifying seed
set.seed(333422)#simulating trajectory
for (i in 2:(njumps+1))
time[i]<- time[i-1]+round((-1/lambda)*log(runif(1)),2)Plotting trajectory
# type="n" draws empty frame with no graph
plot(time, N, type="n", xlab="Time", ylab="State",
panel.first = grid())
segments(time[-length(time)], N[-length(time)],
time[-1]-0.07, N[-length(time)], lwd=2, col="blue")
points(time, N, pch=20, col="blue")
points(time[-1], N[-length(time)], pch=1, col="blue")
Grafik tersebut menunjukkan hasil simulasi proses Poisson dengan
parameter laju λ=2. Grafik berbentuk tangga karena pada proses Poisson
nilai proses akan tetap konstan selama tidak ada kejadian, kemudian naik
satu unit setiap kali terjadi suatu kejadian baru. Terlihat bahwa proses
dimulai dari state 0 pada waktu t=0, lalu meningkat secara bertahap
hingga mencapai state 20 pada sekitar waktu t≈12. Jarak antar lompatan
tidak sama karena waktu kedatangan antar kejadian bersifat acak dan
mengikuti distribusi eksponensial. Pada beberapa interval, seperti
sekitar waktu 6 hingga 9, tidak banyak terjadi kejadian sehingga grafik
tampak datar lebih panjang. Sebaliknya, di sekitar waktu 9 hingga 12
terjadi banyak lompatan dalam waktu singkat yang menunjukkan kejadian
muncul lebih rapat.
3.2. Uniform Order Statistics
#specifying parameters
t<- 10
lambda<- 2#specifying seed
set.seed(32114)#generating N(t)
njumps<- rpois(1,lambda*t)#defining states
N<- 0:njumps#generating N(t) standard uniforms
u<- c()
u[1]<- 0
for(i in 2:(njumps+1))
u[i]<- runif(1)#computing event times
time<- t*sort(u)#plotting trajectory
plot(time, N, type="n", xlab="Time", ylab="State",
panel.first = grid())
segments(time[-length(time)], N[-length(time)],
time[-1]-0.07, N[-length(time)], lwd=2, col="blue")
points(time, N, pch=20, col="blue")
points(time[-1], N[-length(time)], pch=1, col="blue")
Grafik adalah realisasi dari suatu proses Poisson di mana, ketika
diketahui bahwa terdapat sejumlah n kejadian hingga waktu t, maka waktu
terjadinya kejadian-kejadian tersebut tersebar seperti order statistics
dari distribusi Uniform(0, t). Hal ini berarti meskipun secara
keseluruhan proses mengikuti dinamika acak (random arrivals), jika
jumlah total kejadian dikunci, maka posisi waktu kejadian tersebut akan
tampak tersebar acak namun merata di sepanjang interval waktu
pengamatan. Pada grafik, titik-titik kejadian tersebar di sepanjang
sumbu waktu tanpa pola deterministik tertentu, tetapi tetap menunjukkan
akumulasi yang meningkat secara bertahap. Variasi jarak antar kejadian
(interarrival times) yang terlihat tidak seragam konsisten dengan sifat
acak dari distribusi uniform setelah diurutkan, di mana beberapa
kejadian bisa tampak berdekatan sementara yang lain berjauhan.
3.3.
3.1. Occurence of Earthquake
library(readr)## Warning: package 'readr' was built under R version 4.4.3library(dplyr)## Warning: package 'dplyr' was built under R version 4.4.3##
## Attaching package: 'dplyr'## The following objects are masked from 'package:stats':
##
## filter, lag## The following objects are masked from 'package:base':
##
## intersect, setdiff, setequal, uniongetwd()## [1] "C:/INDY"#storing data
cols <- c(
"DATE", "TIME", "ET", "GT", "MAG", "M",
"LAT", "LON", "DEPTH", "Q", "EVID", "NPH", "NGRM"
)
eq.data <- read_table(
"C:/Users/KOKO MEME/Documents/data.txt",
col_names = cols,
skip = 1
)##
## ── Column specification ────────────────────────────────────────────────────────
## cols(
## DATE = col_date(format = ""),
## TIME = col_time(format = ""),
## ET = col_character(),
## GT = col_character(),
## MAG = col_double(),
## M = col_character(),
## LAT = col_double(),
## LON = col_double(),
## DEPTH = col_double(),
## Q = col_character(),
## EVID = col_double(),
## NPH = col_double(),
## NGRM = col_double()
## )head(eq.data)## # A tibble: 6 × 13
## DATE TIME ET GT MAG M LAT LON DEPTH Q EVID
## <date> <time> <chr> <chr> <dbl> <chr> <dbl> <dbl> <dbl> <chr> <dbl>
## 1 2012-01-02 21:11:29 eq l 3 l 31.7 -116. 8.6 C 11049077
## 2 2012-01-03 14:18:56 eq l 4.14 l 33.2 -119. 18 B 11049285
## 3 2012-01-06 02:47:51 eq l 3.08 l 34.0 -116. 4.4 A 11050165
## 4 2012-01-06 06:50:10 eq l 3.25 l 36.0 -118. 5 A 11050245
## 5 2012-01-06 16:58:54 eq l 3.47 l 33.7 -116. 4.5 A 11050541
## 6 2012-01-08 01:59:54 eq l 3.51 l 33.7 -116. 6.4 A 11051197
## # ℹ 2 more variables: NPH <dbl>, NGRM <dbl>#creating date-time variable
datetime<- as.POSIXct(paste(as.Date(eq.data$DATE),
eq.data$TIME))#computing lag
datetime.lag<- c(0,head(datetime, -1))#computing interarrival times (in hours)
int.time<- (as.numeric(datetime)-as.numeric(datetime.lag))/
3600#removing first value
int.time<- int.time[-1]#removing immediate aftershocks (within 3 hours)
int<- int.time[int.time>3]#plotting histogram
hist(int, main="", col="dark magenta", xlab="Interarrival
Time")
Histogram menggambarkan distribusi interarrival times (jarak waktu antar
gempa) setelah data dibersihkan dari kemungkinan aftershocks. Secara
visual, distribusi tersebut tampak miring ke kanan (right-skewed), di
mana sebagian besar waktu antar kejadian bernilai kecil (dekat nol), dan
frekuensinya menurun secara bertahap seiring bertambahnya waktu. Pola
ini merupakan ciri khas dari distribusi eksponensial. Hal ini
menunjukkan bahwa gempa bumi cenderung terjadi dengan interval waktu
yang pendek lebih sering, sementara interval yang panjang relatif jarang
terjadi. Ini konsisten dengan sifat memoryless dari distribusi
eksponensial, di mana peluang terjadinya gempa berikutnya tidak
bergantung pada kapan gempa sebelumnya terjadi. Dengan kata lain,
meskipun beberapa interval panjang tetap muncul (terlihat dari ekor
panjang di kanan histogram), mayoritas kejadian terjadi dalam waktu yang
relatif dekat satu sama lain.
Uji Chi-Square
#binning interarrival times
binned.int<- as.factor(ifelse(int<40,"1",
ifelse(int>=40 & int<80,"2",ifelse(int>=80 & int<120,"3",
ifelse(int>=120 & int<160,"4",ifelse(int>=160 &
int<200,"5",
ifelse(int>=200 & int<240,"6","7")))))))#computing observed frequencies
obs<- table(binned.int)#estimating mean for exponential distribution
mean.est<- mean(int)#computing expected frequencies
exp<- c(1:7)
exp[1]<- length(int)*(1-exp(-40/mean.est))
exp[2]<- length(int)*(exp(-40/mean.est)-exp(-80/mean.est))
exp[3]<- length(int)*(exp(-80/mean.est)-exp(-120/mean.est))
exp[4]<- length(int)*(exp(-120/mean.est)-exp(-160/mean.est))
exp[5]<- length(int)*(exp(-160/mean.est)-exp(-200/mean.est))
exp[6]<- length(int)*(exp(-200/mean.est)-exp(-240/mean.est))
exp[7]<- length(int)*exp(-240/mean.est)
obs## binned.int
## 1 2 3 4 5 6 7
## 339 179 117 49 38 25 42#rounding
round(exp,1)## [1] 318.3 189.9 113.3 67.6 40.3 24.0 35.6#computing chi-squared statistic
chi.sq<- sum((obs-exp)^2/exp)
print(chi.sq)## [1] 8.534049#computing p-value
print(p.value<- 1-pchisq(chi.sq, df=5))## [1] 0.129156Berdasarkan hasil uji chi-squared yang diperoleh, nilai p-value lebih besar dari tingkat signifikansi umum (α = 0.05), maka tidak terdapat cukup bukti untuk menolak hipotesis nol. Dengan demikian, dapat disimpulkan bahwa distribusi interarrival time tidak berbeda secara signifikan dari distribusi eksponensial.
3.2. Scoring Goals
a
lambda_A <- 1/10
lambda_B <- 1/12
lambda_total <- lambda_A + lambda_B
expected_wait <- 1 / lambda_total
expected_wait## [1] 5.454545b
p_A_first <- lambda_A / (lambda_A + lambda_B)
p_B_first <- lambda_B / (lambda_A + lambda_B)
p_A_first## [1] 0.5454545p_B_first## [1] 0.4545455c
# Probability Ties
sum<- 0
for(n in 0:15){
sum<- sum+(30^n)/(factorial(n))^2
}
sum*exp(-11)## [1] 0.1165575# Probability A wins
sum.n<- 0
for (n in 0:15) {
sum.k<-0
for (k in 1:15){
sum.k<- sum.k+6^k/factorial(n+k)
sum.n<- sum.n + 30^n/factorial(n)*sum.k
}
}
sum.n*exp(-11)## [1] 7.077008# Probability B wins
sum.n<- 0
for (n in 0:15) {
sum.k<-0
for (k in 1:15){
sum.k<- sum.k+5^k/factorial(n+k)
sum.n<- sum.n + 30^n/factorial(n)*sum.k
}
}
sum.n*exp(-11)## [1] 4.324333