R：restricted cubic splineの作図について

Question

前提・実現したいこと

初めて質問させていただきます。よろしくお願いいたします。
現在、Rでrmsパッケージを用いてcox回帰分析を行っております。
説明変数が非線形であり、restricted cubic splineを実行しようと考えています。

作図に関しては下記サイトのようにな2群に分けて行いたいのですが、様々なwebサイトを探してみてもやり方が見つかりませんでした。
https://www.scielo.br/scielo.php?pid=S0034-70942020005003213&script=sci_arttext#aff1

もし方法等ご存知な方がおられましたらご教示いただけますと幸いです。よろしくお願いいたします。

該当のソースコード

install.packages("rms")
library(rms)
install.packages("datadist")
ddist1 <- datadist(df1)
options(datadist='ddist1')
k <- with(df1, quantile(df1$heart rate, c(.05, 0.25, .75, .95)))
k

fit <- cph(Surv(time, death) ~ rcs(heart rate, k) + category_Age + category_BP + category_BUN + category_Na + copd, data = df1, nk = 5)
plot(Predict(fit, heart rate), xlab = "heart rate", ylab = "Relative Risk", lty = 1, lwd = 2)

2群に分けずに全集団で作図するところまでは何とかできたのですが、2群に分けて作図する方法が分かりませんでした。

Answer 1

例えば、category_Ageが40歳未満、40歳以上といったような2群のカテゴリとした場合、

plot(Predict(fit, heart rate, category_Age), xlab = "heart rate", ylab = "Relative Risk")

でできた気がします。

さらに

contrast(fit, list(category_Age =0, heart rate=c(50,80,120)),list(category_Age =1, heart rate=c(50,80,120)))

みたいにすれば、Ageカテゴリ（0群と1群）群間において、50、80、120のそれぞれの心拍数における有意差検定ができたような気がします。

間違ってたらすみません。

Answer 2

次のような方針で考えるとよいと思います。

1. 「全集団で作図するところまでは何とかできた」部分を一つの関数にまとめる。
2. tapply(データ, 属性,上で定義した関数)でグループごとの解析を実行する。

参考:

• R-Source
• Rでクラスごとの記述統計を効率的に算出するお話 | Takuro FUJITA@NIT, Fukui College