Multiple regression model

set.seed(123)
n <- 100  
x1 <- rnorm(n, mean = 50, sd = 10)  # Quantitative predictor
x2 <- sample(0:1, n, replace = TRUE)  # Dichotomous predictor
y <- 2*x1 + 0.5*x1^2 + 3*x2 + 0.5*x1*x2 + rnorm(n, mean = 0, sd = 5)  # Dependent variable
data <- data.frame(y, x1, x2)

# Multiple regression model
model <- lm(y ~ x1 + I(x1^2) + x2 + x1:x2, data = data)
summary(model)
## 
## Call:
## lm(formula = y ~ x1 + I(x1^2) + x2 + x1:x2, data = data)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -8.0471 -3.1952 -0.7353  2.2315 14.6579 
## 
## Coefficients:
##               Estimate Std. Error t value Pr(>|t|)    
## (Intercept) -15.306248  10.796547  -1.418    0.160    
## x1            2.720654   0.437101   6.224 1.30e-08 ***
## I(x1^2)       0.491693   0.004397 111.825  < 2e-16 ***
## x2           -0.811403   5.691925  -0.143    0.887    
## x1:x2         0.608250   0.109719   5.544 2.65e-07 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 4.557 on 95 degrees of freedom
## Multiple R-squared:  0.9999, Adjusted R-squared:  0.9999 
## F-statistic: 2.891e+05 on 4 and 95 DF,  p-value: < 2.2e-16

Residual analysis

# Residual analysis
par(mfrow = c(2, 2))
plot(model)

The points in the QQ plot approximately follow a straight line, indicating that the residuals are normally distributed.