Understanding the Johnson-Neyman algorithm for finding the region of significance in a trivariate correlation analysis

Question

I'm pretty new to the whole statistics world, so excuse me for being a novice in advance. I'm trying to understand Andrew F. Hayes's PROCESS macro in SPSS [1]. To my best understanding, the macro uses the Johnson-Neyman [JN] method [2] to partition a quantitative moderator/mediator value into some regions of significance. And then uses the ANOVA method to analyze the bivariate correlation of the two other variables grouped based on the calculated intervals of the moderator value. Something like a non-heuristic classification method.

Here Carden et al, [3] have implemented the JN method in an Excell template called CAHOST (Download). Other implementations of the algorithm should be available in PyProcessMacro, processr, and ProcessR. Now as I have explained here, my goal is to understand the math behind JN and maybe implement a macro for Open/LibreOffice Calc. I would appreciate it if you could walk me through the math behind the JN method, for example as explained here. Consider the first equation:

$$Y_i = \gamma_0 + \gamma_1 X_i + \gamma_2 M_i + \gamma_3 X_i M_i + \epsilon_i \tag{1}$$

to be the assumed trivariate linear casualty-moderation/mediation between $X$ as the focal predictor, $Y$ is the response variable and $M$ is the moderator/mediator. The $\gamma_j$ coefficients should be constant scalars and $\epsilon_i$ is the regression error for the $i^{th}$ sample. They have written the equation in matrix form for the whole sample group as

$$\hat{Y} = \hat{\gamma}_0 + \hat{\gamma}_1 X + \hat{\gamma}_2 M + \hat{\gamma}_3 X M \tag{2}$$

which I do not understand. For example, what are the dimensions of the matrices in the above equation? Shouldn't the $\hat{\gamma}_j$s be constant scalars? and where did the error go?

References:

1. Hayes, A. F. (2013b). Introduction to Mediation, Moderation, and Conditional Process Analysis: A Regression-Based Approach. New York, NY: The Guilford Press.
2. Johnson, P. O., and Neyman, J. (1936). Tests of certain linear hypotheses and their applications to some educational problems. Stat. Res. Mem. 1, 57–39.
3. Carden, S. W., Holtzman, N. S., & Strube, M. J. (2017). CAHOST: An Excel workbook for facilitating the Johnson-Neyman technique for two-way interactions in multiple regression. Frontiers in Psychology, 8, 1293.

Answer 1

As an admitted statistical novice, it seems there are a few concepts from the statistical paradigm not yet grasped, as well as a few errors in interpretation.

Think of equation (1) as the "true" model. The $\gamma_j$ are indeed constant scalars, but we don't know what they are. The error term $\epsilon$ is there to indicate that the responses will not match the deterministic portion of the model exactly. There will be variance.

Neither equations (1) or (2) are written in matrix form. All quantities are still scalar.

The $\hat{\gamma_j}$ terms represent estimates of the "true" but unknown $\gamma_j$ terms. Once you've collected some data and estimated them, then they can be regarded as constant scalars. But thinking more generally, before the data is collected, they are random variables that depend on the data. The values calculated via least squares represent only one possible realization. This perspective is key to understanding the remaining derivation of the JN method, which uses the variance of the $\hat{\gamma_j}$ terms.

Once data is collected and $\hat{\gamma_j}$ are estimated, equation (2) can be used to predict responses for new observations. It is meant to represent the average response for given predictor/moderator values (technically it is a conditional expectation), and the error term is assumed to have zero mean, so it is not included.

I hope this gets you started.