vgam2

VGAM

# パッケージのインストール
#install.packages("VGAM")

# パッケージのロード
library(VGAM)

# データ生成
set.seed(42)
age <- seq(20, 70, length.out = 100)
diastolic_bp <- 70 + 0.5 * age + rnorm(100, mean = 0, sd = 10)
systolic_bp <- 100 + 0.7 * age + rnorm(100, mean = 0, sd = 15)

# データフレームにまとめる
data <- data.frame(age = age, diastolic_bp = diastolic_bp, systolic_bp = systolic_bp)

# データのプロット
plot(data$age, data$diastolic_bp, pch = 1, col = "blue", xlab = "Age", ylab = "Blood Pressure", main = "Generated Blood Pressure Data", ylim = c(50,200))
points(data$age, data$systolic_bp, pch = 4, col = "red")
legend("topleft", legend = c("Diastolic BP (〇)", "Systolic BP (×)"), pch = c(1, 4), col = c("blue", "red"))

# 同時に回帰スプラインを適用
model_vglm <- vglm(cbind(diastolic_bp, systolic_bp) ~ bs(age, df = 4), family = gaussianff(), data = data)

# 予測値の取得
pred_vglm <- predict(model_vglm, newdata = data)

# プロット
plot(data$age, data$diastolic_bp, pch = 1, col = "blue", xlab = "Age", ylab = "Blood Pressure", main = "Regression Splines with vglm",ylim=c(50,250))
points(data$age, data$systolic_bp, pch = 4, col = "red")
lines(data$age, pred_vglm[,1], col = "blue", lwd = 2)
lines(data$age, pred_vglm[,3], col = "red", lwd = 2)
legend("topleft", legend = c("Diastolic BP (〇)", "Systolic BP (×)", "Diastolic BP Spline Fit", "Systolic BP Spline Fit"), pch = c(1, 4, NA, NA), lty = c(NA, NA, 1, 1), col = c("blue", "red", "blue", "red"), lwd = c(NA, NA, 2, 2))

# 同時にベクトルスムージングスプラインを適用
model_vgam_gam <- vgam(cbind(diastolic_bp, systolic_bp) ~ s(age, df = 4), family = gaussianff(), data = data)

# 予測値の取得
pred_vgam_gam <- predict(model_vgam_gam, newdata = data)

# プロット
plot(data$age, data$diastolic_bp, pch = 1, col = "blue", xlab = "Age", ylab = "Blood Pressure", main = "Vector Smoothing Splines with vgam",ylim=c(50,250))
points(data$age, data$systolic_bp, pch = 4, col = "red")
lines(data$age, pred_vgam_gam[,1], col = "blue", lwd = 2)
lines(data$age, pred_vgam_gam[,3], col = "red", lwd = 2)
legend("topleft", legend = c("Diastolic BP (〇)", "Systolic BP (×)", "Diastolic BP GAM Fit", "Systolic BP GAM Fit"), pch = c(1, 4, NA, NA), lty = c(NA, NA, 1, 1), col = c("blue", "red", "blue", "red"), lwd = c(NA, NA, 2, 2))

# 同時にPスプラインを適用
model_vgam_ps <- vgam(cbind(diastolic_bp, systolic_bp) ~ s(age, bs = "ps", df = 4), family = gaussianff(), data = data)

# 予測値の取得
pred_vgam_ps <- predict(model_vgam_ps, newdata = data)

# プロット
plot(data$age, data$diastolic_bp, pch = 1, col = "blue", xlab = "Age", ylab = "Blood Pressure", main = "P-Splines with vgam",ylim=c(50,250))
points(data$age, data$systolic_bp, pch = 4, col = "red")
lines(data$age, pred_vgam_ps[,1], col = "blue", lwd = 2)
lines(data$age, pred_vgam_ps[,3], col = "red", lwd = 2)
legend("topleft", legend = c("Diastolic BP (〇)", "Systolic BP (×)", "Diastolic BP P-Spline Fit", "Systolic BP P-Spline Fit"), pch = c(1, 4, NA, NA), lty = c(NA, NA, 1, 1), col = c("blue", "red", "blue", "red"), lwd = c(NA, NA, 2, 2))