我目前正在研究一个混合模型,其中包括分类变量和度量成本变量的多个交互,以预测度量输出。由于众多的分类变量,我对成本和产出(斜率)之间的线性关系的变化非常感兴趣。
因为我有多个与成本交互的分类协变量(并且还有随机效应),所以我寻找一种简单的方法来计算这些斜率的边际效应。我不是在寻找对输出的预测或期望,而是希望看到“边际斜率”。
我在下面包含了一些示例代码
library(housingData)
library(tidyverse)
library(lme4)
#> Loading required package: Matrix
#>
#> Attaching package: 'Matrix'
#> The following objects are masked from 'package:tidyr':
#>
#> expand, pack, unpack
set.seed(12345)
df <- housing %>%
filter(complete.cases(.)) %>%
mutate(house_type = factor(sample(c("A", "B", "C"), nrow(.), replace = TRUE)),
timeperiod = as.factor(if_else(lubridate::year(time) < 2012,
"before_2012",
"2012_and_after")))
# model categorical effects on slope, after main effect has been fitted
lm1 <- lm(medSoldPriceSqft ~ 0 + medListPriceSqft, data = df)
coef(lm1)
#> medListPriceSqft
#> 0.9137324
df2 <- df %>%
modelr::add_residuals(lm1, "medSoldPriceSqft_after_list_price")
mod2 <- lmer(medSoldPriceSqft_after_list_price ~ 0 + medListPriceSqft:house_type
+ medListPriceSqft:timeperiod +
(0 + medListPriceSqft:house_type + medListPriceSqft:timeperiod| state), data = df2)
fixef(mod2)
#> medListPriceSqft:house_typeA medListPriceSqft:house_typeB
#> 0.02819789 0.01725657
#> medListPriceSqft:house_typeC medListPriceSqft:timeperiodbefore_2012
#> 0.01733967 0.01526881
ranef(mod2)$state
#> medListPriceSqft:house_typeA medListPriceSqft:house_typeB
#> AL -0.0332371439 -0.0342472679
#> AR -0.0161834913 -0.0197485094
#> AZ -0.0455633279 -0.0540328340
#> CA -0.0073736100 -0.0069451237
#> CO -0.0998324113 -0.2144934915
#> CT -0.0068844784 -0.0064218278
#> DC 0.1323688987 0.1625294726
#> DE 0.0194554114 0.0240963193
#> FL -0.0937513315 -0.1128399131
#> GA -0.0404177163 -0.0447775564
#> IA -0.0018353661 -0.0065783305
#> ID 0.0715311714 0.0827129047
#> IL 0.0364040836 0.0519384804
#> IN 0.3815562061 0.5127737560
#> KS 0.3195998033 0.3735174115
#> KY 0.0180534996 0.0284589475
#> LA -0.0034512710 0.0002398266
#> MA 0.0350797588 0.0495612406
#> MD -0.0061851180 0.0030358059
#> ME -0.1205435796 -0.1434406842
#> MI -0.0104773953 -0.0097225674
#> MN 0.0658291099 0.0815773829
#> MO 0.1658471055 0.2092785032
#> MS 0.0757674725 0.0957259331
#> MT -0.2226956530 -0.2848157760
#> NC -0.0801849629 -0.0930100202
#> ND 0.1471279546 0.1632516890
#> NE -0.0067148446 -0.0091668294
#> NH -0.0288545182 -0.0318590020
#> NJ 0.0419945037 0.0452393238
#> NM -0.3402646109 -0.4196645043
#> NV -0.0319233561 -0.0379518075
#> NY -0.0458677598 -0.0460468454
#> OH 0.0053817946 0.0082642996
#> OK -0.0197739286 -0.0218087721
#> OR -0.0312838479 -0.0325376041
#> PA 0.0548995292 0.0837696293
#> RI -0.0301300132 -0.0360929264
#> SC -0.0627597227 -0.0749119310
#> TN 0.0327156573 0.0397181302
#> TX 0.0813191504 0.0826967005
#> UT -0.0078590502 -0.0157365323
#> VA -0.0007437153 0.0095222182
#> VT -0.1811751876 -0.2204967507
#> WA -0.0090309685 -0.0112995414
#> WI -0.0244199376 -0.0283874222
#> WV -0.0755127962 -0.0908736075
#> medListPriceSqft:house_typeC medListPriceSqft:timeperiodbefore_2012
#> AL -0.0161539511 -0.033273312
#> AR -0.0143981558 -0.002542357
#> AZ -0.0354905346 -0.022865789
#> CA -0.0023188159 -0.016206762
#> CO -0.1891388131 -0.044702463
#> CT -0.0025396683 -0.007847727
#> DC 0.0915280219 0.147291653
#> DE 0.0138827021 0.019724991
#> FL -0.0701567701 -0.084079270
#> GA -0.0260708779 -0.036974280
#> IA -0.0102049733 0.015404515
#> ID 0.0420144801 0.079569678
#> IL 0.0595382970 -0.072987720
#> IN 0.4171418424 -0.035008222
#> KS 0.1945214096 0.349026709
#> KY 0.0298359996 -0.024003880
#> LA 0.0102208245 -0.036621369
#> MA 0.0494285509 -0.034373164
#> MD 0.0074462978 -0.007041631
#> ME -0.0886152719 -0.076925785
#> MI -0.0009559257 -0.024871313
#> MN 0.0334029227 0.132315832
#> MO 0.1447632248 0.074557543
#> MS 0.0624490008 0.054095381
#> MT -0.2229587595 0.002527347
#> NC -0.0581766037 -0.062570993
#> ND 0.0566945904 0.254531126
#> NE -0.0065563090 -0.003424448
#> NH -0.0121567614 -0.046952309
#> NJ 0.0457664635 -0.067748718
#> NM -0.2712390810 -0.211903077
#> NV -0.0242007544 -0.018678370
#> NY -0.0430582366 0.108457835
#> OH 0.0151824401 -0.038200489
#> OK -0.0100685546 -0.023278317
#> OR -0.0151429366 -0.036220790
#> PA 0.0781025694 -0.026947879
#> RI -0.0300959629 0.012527652
#> SC -0.0372963600 -0.092552679
#> TN 0.0308868734 -0.001932841
#> TX 0.0381596134 0.076102496
#> UT -0.0305654079 0.070222766
#> VA 0.0037802116 0.046899465
#> VT -0.1400341335 -0.118380902
#> WA -0.0022661974 -0.027577840
#> WI -0.0106677175 -0.046255988
#> WV -0.0542188035 -0.060304313
# make marginal effects accessible
## Average increase in medListPriceSqft-Slope for house_typeA and timeperiodbefore_2012
fixef(mod2)["medListPriceSqft:house_typeA"] + fixef(mod2)["medListPriceSqft:timeperiodbefore_2012"]
#> medListPriceSqft:house_typeA
#> 0.0434667
## Average increase in medListPriceSqft-Slope for house_typeA ()
prop_before_2012 <- table(df2$timeperiod)["before_2012"] / nrow(df2)
fixef(mod2)["medListPriceSqft:house_typeA"] +
fixef(mod2)["medListPriceSqft:timeperiodbefore_2012"] * prop_before_2012
#> medListPriceSqft:house_typeA
#> 0.03564273
## Average increase in medListPriceSqft-Slope for timeperiodbefore_2012
# ...
## Average increase in medListPriceSqft-Slope for house_typeA () for state == "WA"
# ...
如您所见,这很快就会变得乏味,我怎样才能以编程方式做到这一点?包装{margins}
或{ggeffects}
似乎仅涵盖输出的边际效应?