library(rockchalk)

R.rc<-lm(Y ~ X*Z, data=Moderation)

summary(R.rc)
## 
## Call:
## lm(formula = Y ~ X * Z, data = Moderation)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -1.68679 -0.39206  0.01796  0.42799  2.19791 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept) -0.12376    0.55564  -0.223 0.823970    
## X            0.46893    0.12958   3.619 0.000377 ***
## Z            0.09560    0.13146   0.727 0.467960    
## X:Z          0.07861    0.03005   2.616 0.009592 ** 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.6799 on 196 degrees of freedom
## Multiple R-squared:  0.7546, Adjusted R-squared:  0.7508 
## F-statistic: 200.9 on 3 and 196 DF,  p-value: < 2.2e-16

Plots

PS.rc<-plotSlopes(R.rc, plotx="X", modx="Z", modxVals="std.dev")

TS.rc<-testSlopes(PS.rc)
## Values of Z OUTSIDE this interval:
##         lo         hi 
## -36.879489  -1.571709 
## cause the slope of (b1 + b2*Z)X to be statistically significant
plot(TS.rc)

jtools plot

library(jtools)
## 
## Attaching package: 'jtools'
## The following object is masked from 'package:rockchalk':
## 
##     standardize
MOD3<-lm(Y ~ X*Z, data=Moderation)

johnson_neyman(model = MOD3, pred = X, modx = Z)
## JOHNSON-NEYMAN INTERVAL 
## 
## When Z is OUTSIDE the interval [-36.88, -1.57], the slope of X is p <
## .05.
## 
## Note: The range of observed values of Z is [1.00, 7.00]