傾向スコア
BABA Yoshihiko
2021/10/4
論文 Masaki S and Kawamoto T (2019)
Masaki S & Kawamoto T (2019) Comparison of long-term outcomes between enteral nutrition via gastrostomy and total parenteral nutrition in older persons with dysphagia: A propensity-matched cohort study. PLoS One, 14(10), e0217120.
library(readxl)
dfPropensityScore <- read_excel("Dataset_PEG_TPN_Dryad.xlsx")
dfPropensityScore$PEG <- factor(dfPropensityScore$PEG,
levels = c(1,0),
labels = c("PEG","TPN"))
dfPropensityScore$PORT <- as.factor(dfPropensityScore$PORT)
dfPropensityScore$NTCVC <- as.factor(dfPropensityScore$`NT-CVC`)
dfPropensityScore$PICC <- as.factor(dfPropensityScore$PICC)
dfPropensityScore$sex <- as.factor(dfPropensityScore$sex)
dfPropensityScore$CI <- as.factor(dfPropensityScore$CI)
dfPropensityScore$dement <- as.factor(dfPropensityScore$dement)
dfPropensityScore$NMD <- as.factor(dfPropensityScore$NMD)
dfPropensityScore$asp <- as.factor(dfPropensityScore$asp)
dfPropensityScore$IHD <- as.factor(dfPropensityScore$IHD)
dfPropensityScore$lung <- as.factor(dfPropensityScore$lung)
dfPropensityScore$liver <- as.factor(dfPropensityScore$liver)
dfPropensityScore$CKD <- as.factor(dfPropensityScore$CKD)表１の作成
library(tableone)
CreateTableOne(data = dfPropensityScore,
strata="PEG",
vars=c("PORT","NTCVC","PICC", "age", "sex","CI","dement","NMD","asp","IHD","CHF","lung","liver","CKD","Alb", "TLC", "TC", "Hb", "CRP"),
factorVars=c("PORT","NTCVC","PICC", "sex","CI","dement","NMD","asp","IHD","CHF","lung","liver","CKD"))## Stratified by PEG
## PEG TPN p test
## n 180 73
## PORT = 1 (%) 0 ( 0.0) 28 (38.4) <0.001
## NTCVC = 1 (%) 0 ( 0.0) 26 (35.6) <0.001
## PICC = 1 (%) 0 ( 0.0) 19 (26.0) <0.001
## age (mean (SD)) 81.56 (9.84) 86.77 (6.72) <0.001
## sex = 1 (%) 109 (60.6) 45 (61.6) 0.985
## CI = 1 (%) 107 (59.4) 26 (35.6) 0.001
## dement = 1 (%) 57 (31.7) 45 (61.6) <0.001
## NMD = 1 (%) 10 ( 5.6) 4 ( 5.5) 1.000
## asp = 1 (%) 73 (40.6) 21 (28.8) 0.106
## IHD = 1 (%) 31 (17.2) 16 (21.9) 0.489
## CHF = 1 (%) 70 (38.9) 37 (50.7) 0.114
## lung = 1 (%) 12 ( 6.7) 7 ( 9.6) 0.592
## liver = 1 (%) 9 ( 5.0) 6 ( 8.2) 0.491
## CKD = 1 (%) 29 (16.1) 24 (32.9) 0.005
## Alb (mean (SD)) 3.24 (0.57) 2.82 (0.53) <0.001
## TLC (mean (SD)) 1383.84 (713.18) 1203.11 (682.95) 0.073
## TC (mean (SD)) 160.19 (38.06) 145.58 (43.21) 0.009
## Hb (mean (SD)) 11.29 (1.90) 10.21 (2.10) <0.001
## CRP (mean (SD)) 2.03 (2.77) 2.92 (3.04) 0.026MathcIt
まず、欠測値があるため、論文では多重代入法で代入していましたが、ここでは欠測値があるデータを削除します。
#library(mice)
#micePS <- mice(dfPropensityScore)
#dfPropensityScore <- complete(micePS, 1)
dfPropensityScore <- na.omit(dfPropensityScore)次に、傾向スコアマッチングを実行します。
library("MatchIt")
miPropensityScore <- matchit(
PEG ~ age + sex + CI + dement + NMD + asp + IHD + CHF + lung + liver + CKD + Alb + TLC + TC + Hb + CRP,
data = dfPropensityScore,
method = "nearest",
distance = "glm",
caliper = 0.05,
ratio = 1)
summary(miPropensityScore)##
## Call:
## matchit(formula = PEG ~ age + sex + CI + dement + NMD + asp +
## IHD + CHF + lung + liver + CKD + Alb + TLC + TC + Hb + CRP,
## data = dfPropensityScore, method = "nearest", distance = "glm",
## caliper = 0.05, ratio = 1)
##
## Summary of Balance for All Data:
## Means Treated Means Control Std. Mean Diff. Var. Ratio eCDF Mean
## distance 0.4538 0.2081 1.1758 1.2303 0.3177
## age 86.8750 81.1488 0.9077 0.4120 0.1248
## sex0 0.3594 0.4048 -0.0946 . 0.0454
## sex1 0.6406 0.5952 0.0946 . 0.0454
## CI0 0.6250 0.3988 0.4672 . 0.2262
## CI1 0.3750 0.6012 -0.4672 . 0.2262
## dement0 0.3750 0.6845 -0.6393 . 0.3095
## dement1 0.6250 0.3155 0.6393 . 0.3095
## NMD0 0.9375 0.9524 -0.0615 . 0.0149
## NMD1 0.0625 0.0476 0.0615 . 0.0149
## asp0 0.7188 0.5893 0.2879 . 0.1295
## asp1 0.2812 0.4107 -0.2879 . 0.1295
## IHD0 0.7969 0.8274 -0.0758 . 0.0305
## IHD1 0.2031 0.1726 0.0758 . 0.0305
## CHF 0.5000 0.3869 0.2262 . 0.1131
## lung0 0.9219 0.9286 -0.0250 . 0.0067
## lung1 0.0781 0.0714 0.0250 . 0.0067
## liver0 0.9062 0.9464 -0.1378 . 0.0402
## liver1 0.0938 0.0536 0.1378 . 0.0402
## CKD0 0.6562 0.8393 -0.3854 . 0.1830
## CKD1 0.3438 0.1607 0.3854 . 0.1830
## Alb 2.8359 3.2304 -0.7440 0.8398 0.1300
## TLC 1202.6469 1364.2196 -0.2314 0.9319 0.0977
## TC 146.5625 160.1012 -0.3073 1.3330 0.0942
## Hb 10.1609 11.2750 -0.5118 1.3217 0.1277
## CRP 2.6316 2.1298 0.1876 0.8906 0.1114
## eCDF Max
## distance 0.4993
## age 0.2857
## sex0 0.0454
## sex1 0.0454
## CI0 0.2262
## CI1 0.2262
## dement0 0.3095
## dement1 0.3095
## NMD0 0.0149
## NMD1 0.0149
## asp0 0.1295
## asp1 0.1295
## IHD0 0.0305
## IHD1 0.0305
## CHF 0.1131
## lung0 0.0067
## lung1 0.0067
## liver0 0.0402
## liver1 0.0402
## CKD0 0.1830
## CKD1 0.1830
## Alb 0.3058
## TLC 0.2202
## TC 0.2113
## Hb 0.3006
## CRP 0.2805
##
##
## Summary of Balance for Matched Data:
## Means Treated Means Control Std. Mean Diff. Var. Ratio eCDF Mean
## distance 0.3929 0.3907 0.0105 1.0127 0.0062
## age 86.3878 86.3878 0.0000 1.0275 0.0133
## sex0 0.3673 0.3673 0.0000 . 0.0000
## sex1 0.6327 0.6327 0.0000 . 0.0000
## CI0 0.6327 0.6122 0.0422 . 0.0204
## CI1 0.3673 0.3878 -0.0422 . 0.0204
## dement0 0.3878 0.4490 -0.1265 . 0.0612
## dement1 0.6122 0.5510 0.1265 . 0.0612
## NMD0 0.9388 0.9184 0.0843 . 0.0204
## NMD1 0.0612 0.0816 -0.0843 . 0.0204
## asp0 0.6531 0.6939 -0.0908 . 0.0408
## asp1 0.3469 0.3061 0.0908 . 0.0408
## IHD0 0.7959 0.8367 -0.1015 . 0.0408
## IHD1 0.2041 0.1633 0.1015 . 0.0408
## CHF 0.4694 0.4898 -0.0408 . 0.0204
## lung0 0.9592 0.8980 0.2281 . 0.0612
## lung1 0.0408 0.1020 -0.2281 . 0.0612
## liver0 0.9592 0.9184 0.1400 . 0.0408
## liver1 0.0408 0.0816 -0.1400 . 0.0408
## CKD0 0.7347 0.7347 0.0000 . 0.0000
## CKD1 0.2653 0.2653 0.0000 . 0.0000
## Alb 2.9286 2.9776 -0.0924 0.9234 0.0299
## TLC 1168.4143 1261.6408 -0.1335 1.1308 0.0743
## TC 145.5918 151.1429 -0.1260 1.8412 0.0597
## Hb 10.2122 10.2551 -0.0197 2.0379 0.0673
## CRP 2.7278 1.9763 0.2809 1.3598 0.1004
## eCDF Max Std. Pair Dist.
## distance 0.0408 0.0188
## age 0.0612 6.1633
## sex0 0.0000 0.5714
## sex1 0.0000 0.5714
## CI0 0.0204 1.2225
## CI1 0.0204 1.2225
## dement0 0.0612 0.8009
## dement1 0.0612 0.8009
## NMD0 0.0204 0.4215
## NMD1 0.0204 0.4215
## asp0 0.0408 0.9078
## asp1 0.0408 0.9078
## IHD0 0.0408 0.7102
## IHD1 0.0408 0.7102
## CHF 0.0204 0.9388
## lung0 0.0612 0.5323
## lung1 0.0612 0.5323
## liver0 0.0408 0.2801
## liver1 0.0408 0.2801
## CKD0 0.0000 0.3673
## CKD1 0.0000 0.3673
## Alb 0.1224 0.8315
## TLC 0.2653 0.9737
## TC 0.1837 0.9136
## Hb 0.1429 0.9891
## CRP 0.2041 1.0273
##
## Percent Balance Improvement:
## Std. Mean Diff. Var. Ratio eCDF Mean eCDF Max
## distance 99.1 93.9 98.1 91.8
## age 100.0 96.9 89.4 78.6
## sex0 100.0 . 100.0 100.0
## sex1 100.0 . 100.0 100.0
## CI0 91.0 . 91.0 91.0
## CI1 91.0 . 91.0 91.0
## dement0 80.2 . 80.2 80.2
## dement1 80.2 . 80.2 80.2
## NMD0 -37.1 . -37.1 -37.1
## NMD1 -37.1 . -37.1 -37.1
## asp0 68.5 . 68.5 68.5
## asp1 68.5 . 68.5 68.5
## IHD0 -33.8 . -33.8 -33.8
## IHD1 -33.8 . -33.8 -33.8
## CHF 82.0 . 82.0 82.0
## lung0 -814.3 . -814.3 -814.3
## lung1 -814.3 . -814.3 -814.3
## liver0 -1.6 . -1.6 -1.6
## liver1 -1.6 . -1.6 -1.6
## CKD0 100.0 . 100.0 100.0
## CKD1 100.0 . 100.0 100.0
## Alb 87.6 54.4 77.0 60.0
## TLC 42.3 -74.4 24.0 -20.5
## TC 59.0 -112.4 36.6 13.1
## Hb 96.2 -155.3 47.3 52.5
## CRP -49.8 -165.3 9.9 27.2
##
## Sample Sizes:
## Control Treated
## All 168 64
## Matched 49 49
## Unmatched 119 15
## Discarded 0 0関数
- matchit: ライブラリ MatchiIt の関数
引数
- c caliper = 0.aliper = 0.
dfPSMatched <- match.data(miPropensityScore)
library(survival)
library(survminer)## 要求されたパッケージ ggplot2 をロード中です## 要求されたパッケージ ggpubr をロード中ですobjSurv <- survfit( Surv(survival,status) ~ PEG,
data = dfPSMatched)
ggsurvplot(objSurv,
data = dfPSMatched,
risk.table = TRUE,
pval = TRUE,
conf.int = TRUE,
risk.table.y.text.col = TRUE)
dfPSMatched$oral <- as.factor(dfPSMatched$oral)
dfPSMatched$home <- as.factor(dfPSMatched$home)
dfPSMatched$pneumo <- as.factor(dfPSMatched$pneumo)
dfPSMatched$sepsis <- as.factor(dfPSMatched$sepsis)
CreateTableOne(data = dfPSMatched,
strata="PEG",
vars=c("oral", "home", "pneumo", "sepsis"),
factorVars=c("oral", "home", "pneumo", "sepsis"))## Stratified by PEG
## PEG TPN p test
## n 49 49
## oral = 1 (%) 3 ( 6.1) 2 ( 4.1) 1.000
## home = 1 (%) 5 (10.2) 5 (10.2) 1.000
## pneumo = 1 (%) 23 (46.9) 14 (28.6) 0.096
## sepsis = 1 (%) 7 (14.3) 17 (34.7) 0.035Fig 2
#coxPS <- coxph(Surv(survival,status) ~ PEG + sex + strata(PEG),
# data = dfPSMatched)
#
#ggforest(coxPS)
coxPSsex0 <- coxph(Surv(survival,status) ~ PEG,
data = subset(dfPSMatched, sex == 0))
coxPSsex1 <- coxph(Surv(survival,status) ~ PEG,
data = subset(dfPSMatched, sex == 1))
ggforest(coxPSsex0)## Warning in .get_data(model, data = data): The `data` argument is not provided.
## Data will be extracted from model fit.
ggforest(coxPSsex1)## Warning in .get_data(model, data = data): The `data` argument is not provided.
## Data will be extracted from model fit.