View KFAS_crash.R
require(KFAS) # 1.2.9 | |
require(dplyr) | |
require(tidyr) | |
require(ggplot2) | |
require(data.table) | |
require(zoo) | |
##### read dataset ##### | |
# https://catalog.data.gov/dataset/allegheny-county-crash-data | |
# data description |
View KFAS_ARMA.R
require(KFAS) # 1.2.9 | |
require(dplyr) | |
require(tidyr) | |
require(ggplot2) | |
# --- ARIMA(2, 1) with linear trend --- | |
# generate a dataset | |
set.seed(42) | |
t <- 100 | |
y <- arima.sim(n = t, model = list(ar=c(.3, -0.1), ma=.2), sd=.1) + seq(from=1, to=10, length.out = t) |
View bsts.R
require(bsts) # 0.7.1 | |
data(iclaims) # bring the initial.claims data into scope | |
# --- model 1 ---- | |
ss <- AddLocalLinearTrend(list(), initial.claims$iclaimsNSA) | |
ss <- AddSeasonal(ss, initial.claims$iclaimsNSA, nseasons = 52) | |
model1 <- bsts(initial.claims$iclaimsNSA, | |
state.specification = ss, | |
niter = 1000) |
View bsts.R
source("common.R", encoding = "utf-8") | |
df$RP <- calc_RP(df$RP095, df$AP, .95) | |
z <- calc_Z(RP = df$RP, p = df$AP) | |
df <- mutate(df, z1=z$z1, z2=z$z2) | |
df <- mutate(df, z1E=z1*end, z2E=z2*end) | |
ss <- AddLocalLevel(list(), y = df$logPI) # c | |
ss <- AddAr(ss, lags=2, y = df$logPI) # AR(2) | |
# time-varying regression |
View kfas.R
source("common.R", encoding = "utf-8") | |
# ----- KFAS ------ | |
df$RP <- calc_RP(df$RP095, df$AP, .95) | |
z <- calc_Z(RP = df$RP, p = df$AP) | |
df <- mutate(df, z1=z$z1, z2=z$z2) | |
df <- mutate(df, z1E=z1*end, z2E=z2*end) | |
# specify model | |
model3KFAS <- SSModel(logPI ~ SSMtrend(1, Q=NA) + |
View dlm.R
source("common.R", encoding = "utf-8") | |
# ---- dlm ------ | |
# model 3 | |
res <- data.frame() | |
for( a in seq(from=.1, to=.95, by=.05) ){ | |
RP <- calc_RP(RP = df$RP095, p = df$AP, a = a) | |
z <- calc_Z(RP = RP, p = df$AP) | |
z1 <- z$z1 |
View common.R
# ------ common part ---- | |
require(ggplot2) | |
require(dplyr) | |
require(tidyr) | |
require(dlm) # 1.1-4 | |
require(KFAS) # 1.2.9 | |
require(bsts) # 0.7.1 | |
# calculate the reference price |
View mlogit.R
require(mlogit) # 0.2-4 | |
as.mldata <- function(data){ | |
# convert HC dataset | |
# The alternatives are | |
# 1. Gas central heat with cooling (gcc) | |
# 2. Electric central resistence heat with cooling (ecc) | |
# 3. Electric room resistence heat with cooling (erc) | |
# 4. Electric heat pump, which provides cooling also (hpc) | |
# 5. Gas central heat without cooling (gc) |
View Kalman.R
# multivariate normal Kalman filter | |
require(dplyr) | |
require(tidyr) | |
require(ggplot2) | |
require(animation) | |
# ARIMA(1,1) + 線形トレンド の乱数生成 | |
N <- 50 | |
phi1 <- .5 | |
theta1 <- .2 |
View rfm_mcmc_exec.R
library(dplyr) | |
library(tidyr) | |
library(ggplot2) | |
library(rstan) | |
library(loo) | |
library(ggmcmc) | |
# read datasets | |
df <- read.csv("rfm.csv", stringsAsFactors = F) | |
colnames(df)[1] <- "ID" |