What you want is to check the balancing you achieved with PS matching, i.e., whether the descriptives of the treated and the control have reasonably converged.
Matching::MatchBalance provides a rich set of comparative statistics, which helps to draw the respective conclusions.
Here an example using the lalonde data. Note, that you shouldn't apply one-to-one matching when using propensity scores, i.e. the number of matches should not be one, use e.g. M=5. Also a caliper= of e.g. 0.25 standard deviations of maximal acceptable distance for matches should be used to take some account of the common support in both groups.
library(Matching)
data(lalonde)
## Estimate propensity score
fit <- glm(treat ~ age + I(age^2) + educ + I(educ^2) + black +
hisp + married + nodegr + re74 + I(re74^2) + re75 + I(re75^2) +
u74 + u75, family=binomial, data=lalonde)
lalonde$ps <- fit$fitted.values ## store ps in data set
## matching
rr <- with(lalonde, Match(Y=re78, Tr=treat, X=ps, M=5, caliper=.25))
To compare the groups, define a set of comparison variables to easily subset later.
## define comparison variables
(comp <- c('ps', all.vars(fit$call)[1:11]))
# [1] "ps" "treat" "age" "educ" "black" "hisp" "married"
# [8] "nodegr" "re74" "re75" "u74" "u75"
To just compare the groups without any matching, simply aggregate the mean or something else.
aggregate(. ~ treat, lalonde[comp], mean)
# treat ps age educ black hisp married nodegr re74 re75 u74 u75
# 1 0 0.3936490 25.05385 10.08846 0.8269231 0.10769231 0.1538462 0.8346154 2107.027 1266.909 0.7500000 0.6846154
# 2 1 0.4467636 25.81622 10.34595 0.8432432 0.05945946 0.1891892 0.7081081 2095.574 1532.056 0.7081081 0.6000000
To compare the groups before and after matching use Matching::MatchBalance. You may easily use reformulate in the function's formula interface. For the before-after comparison, provide your result rr in match.out=. (To omit console output, use print.level=0.)
bal <- MatchBalance(formul=reformulate(comp[-2], 'treat'), data=lalonde,
match.out=rr)
# ***** (V1) ps *****
# Before Matching After Matching
# mean treatment........ 0.44676 0.43271
# mean control.......... 0.39365 0.4327
# std mean diff......... 45.329 0.010439
#
# mean raw eQQ diff..... 0.054029 0.0017398
# med raw eQQ diff..... 0.060548 0.00085473
# max raw eQQ diff..... 0.10393 0.024826
#
# mean eCDF diff........ 0.13402 0.0043977
# med eCDF diff........ 0.15743 0.0032538
# max eCDF diff........ 0.22443 0.029284
#
# var ratio (Tr/Co)..... 1.2101 0.99967
# T-test p-value........ 1.4842e-06 0.97518
# KS Bootstrap p-value.. < 2.22e-16 0.782
# KS Naive p-value...... 3.7341e-05 0.82407
# KS Statistic.......... 0.22443 0.029284
#
# [...]
#
# Before Matching Minimum p.value: < 2.22e-16
# Variable Name(s): ps Number(s): 1
#
# After Matching Minimum p.value: 0.02
# Variable Name(s): re75 Number(s): 9
Basically, for each variable the function compares the treatment groups before and after matching; mean and standardized mean differences are provided as well as p-values of t-tests and KS tests. At the very bottom, the output states the worst variables according to the t-test p-values.
However, rather than t-tests, in the literature commonly provided are the standardized absolute differences, which in the ideal case should of course be close to zero and may be extracted from the balancing like so:
(sad <- cbind(sapply(bal$BeforeMatching, '[[', 'sdiff.pooled'),
sapply(bal$AfterMatching, '[[', 'sdiff.pooled')) |>
abs() |>
`dimnames<-`(list(comp[-2], c('before', 'after'))))
# ps 47.4345458 0.0104392
# age 10.7277121 7.1516980
# educ 14.1219821 5.5883637
# black 4.3886611 0.0177414
# hisp 17.4561071 0.4161040
# married 9.3640701 2.6733063
# nodegr 30.3986439 3.4347952
# re74 0.2159921 6.2271032
# re75 8.3863254 4.7320098
# u74 9.4140477 6.5848971
# u75 17.6809436 8.5502059
They may easily be plotted using matplot.
matplot(sad, pch=20, xaxt='n', main='Matching success',
ylab='Std. abs. diff.')
arrows(seq_along(comp[-2]), sad[, 1], seq_along(comp[-2]), sad[, 2],
length=.1, col=8)
axis(1, seq_along(comp[-2]), labels=comp[-2])
legend('topright', pch=20, col=1:2, leg=c('before matching', 'after matching'))
Note: R >= 4.1 used.
