# Partial residual plot with confidence intervals in R for multivariate generalized linear regression

I have fit the following regression model:

``````mod <- betareg(connectance ~ fc * size, data = net.land.3000, link = "logit")

summary(mod)

Call:
betareg(formula = connectance ~ fc * size, data = net.land.3000, link = "logit")

Standardized weighted residuals 2:
Min      1Q  Median      3Q     Max
-2.4777 -0.1919 -0.0185  0.1435  2.6781

Coefficients (mean model with logit link):
Estimate Std. Error z value Pr(>|z|)
(Intercept) -2.160853   0.939846  -2.299   0.0215 *
fc           0.062797   0.027081   2.319   0.0204 *
size         0.109050   0.101085   1.079   0.2807
fc:size     -0.005822   0.003427  -1.699   0.0894 .

Phi coefficients (precision model with identity link):
Estimate Std. Error z value Pr(>|z|)
(phi)    85.06      42.32    2.01   0.0445 *

Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Type of estimator: ML (maximum likelihood)
Log-likelihood: 12.88 on 5 Df
Pseudo R-squared: 0.6826
Number of iterations: 35 (BFGS) + 3 (Fisher scoring)
``````

I am trying to plot the relationship between `fc` and `connectance`, with the effects of `size` and the `fc * size` interaction taken into account, and with a 95% confidence interval. From what I understand, this is a partial residual (or marginal effects) plot.

How can I create one in R given that I have multiple independent variables, an interaction, and a non-linear regression? I know the standard method would be to regress the residuals of (y ~ my other independent variables) vs. (my variable of interest ~ other independent variables), but I'm not sure how to handle that with the interaction and non-linear relationship.

I have given it a shot using the `plot_model()` function and got the following:

``````plot_model(mod, type = "pred", terms = "fc")
``````

However, I was expecting the regression line in my partial residual plot to have a slope equal to the model estimate for `fc` (0.06) and an intercept equal to the intercept in my model output (-2.16), which it does not.

To summarize my questions:

1. Is what I want a partial residual plot?
2. If yes, how is that related to the slope of my variable of interest in my model output, if at all?
3. Is `plot_model` giving me the correct results, or is there another function I can use to calculate a regression line and confidence interval?

Sorry if this is very confusing, I am something of a beginner at stats and quite confused. I would appreciate any help, and be happy to provide more information as it is useful!

For reference, here are my data:

``````   connectance       fc size
1    0.3333333 37.96319    8
3    0.2500000 11.33780    8
5    0.3809524 18.16915   13
6    0.5000000 47.88571    5
8    0.2500000 14.02959   10
9    0.1904762 17.87691   13
11   0.2777778 19.11214    9
12   0.2291667 29.03701   14
``````