R语言 PSM之前之后使用Table1输出组间信息对比
2023-06-07 本文已影响0人
王叽叽的小心情
刚开始使用table1,后面使用的是tableone,先备忘录一下
library(tableone)
library(Matching)
library(gtools) # For the logit function
library(lmtest)
library(sandwich)
library(MatchIt)
library(dplyr)
library(table1)
library(printr, quietly=TRUE) # to save table1 plot as csv
library(kableExtra)
library(magick)
pvalue <- function(x, ...) {
# Construct vectors of data y, and groups (strata) g
y <- unlist(x)
g <- factor(rep(1:length(x), times=sapply(x, length)))
if (is.numeric(y)) {
# For numeric variables, perform a standard 2-sample t-test
p <- t.test(y ~ g)$p.value
} else {
# For categorical variables, perform a chi-squared test of independence
p <- chisq.test(table(y, g))$p.value
}
# Format the p-value, using an HTML entity for the less-than sign.
# The initial empty string places the output on the line below the variable label.
c("", sub("<", "<", format.pval(p, digits=3, eps=0.001)))
}
# read data
df <- read.csv("E:\\panel_data.csv", header = 1)
df <- na.omit(df)
# generate the treatment variable
df <- df %>% mutate(trent = ifelse(df$nrty > 2014, 1, 0))
fig_path = "E:\\DID20230606\\"
# before PSM
t_before <- table1(~ df[['ntl']] + poi_comm + poi_trans + lucc_forest + lucc_grass +
poi_publi | treatment, data=df, overall=F,
extra.col=list(`P-value`=pvalue))
# save the figure
t_before[1] %>% save_kable(file = paste(fig_path, "ntl", "_table_before_psm.png"))
# save the numbers
write.csv(as.data.frame(t_before), paste(fig_path, "ntl", "_table_before_psm.csv"),
row.names=TRUE)
# PSM
set.seed(50)
psm <- matchit(data = df,
formula = treatment ~ poi_comm + poi_trans + lucc_forest +
lucc_grass + poi_publi,
method = "nearest",
distance = "logit",
replace = FALSE,
caliper = 0.05)
# after matched
df_matched <- match.data(psm)
t_after <- table1( ~ ntl + poi_comm + poi_trans + lucc_forest + lucc_grass +
poi_publi | treatment, data=df_matched,
overall=F, extra.col=list(`P-value`=pvalue))
# save the figure
t_after[1] %>% save_kable(file = paste(fig_path, "ntl", "_table_after_psm.png"))
# save the numbers
write.csv(as.data.frame(t_after), paste(fig_path, "ntl", "_table_after_psm.csv"),
row.names=TRUE)