# How to calculate marginal effects with multinomial regression for survey data in R

I am trying to calculate average marginal effects for a multinomial logistic regression fitted using the svrepmisc package in R. The margins() function provides solution for binary logistic regression but it does not work for multinomial logistic regression. Is there a package to calculate marginal effects for survey multinomial regression in R? I have included an example illustrating what I have tried to do.

Using margins for svyglm() as follows,

# Binary logistic regression

data(api)
library(survey)
library(margins)

dclus1 <- svydesign(id = ~dnum, weights = ~pw,data = apiclus1, fpc = ~fpc) # survey design
model1 <- svyglm(api00 ~ api99 + sch.wide, design = dclus1) # fit model
margins(model1, design = dclus1) # margins


produces the following output:

Average marginal effects (survey-weighted)
svyglm(formula = api00 ~ api99 + sch.wide, design = dclus1)
api99 sch.wideYes
0.9071       47.48


I tried to use the margins() function for svymultinom() as follows:

# Multinomial logistic regression

library(svrepmisc)

rclus1 <- as.svrepdesign(dclus1) # survey design
model2 <- svymultinom(api00 ~ api99 + sch.wide, design = rclus1)
margins(model2, design = rclus1)


It produces the following error:

Error in model[["call"]] : subscript out of bounds


I am the author of the svrepmisc package that implements svymultinom. I admit that one of the main reasons this package is not on CRAN is that I never went to the trouble to make sure it integrated nicely with other packages like margins, and I also cannot 100% vouch for the accuracy of the variance estimators. Please use with caution.