Tilt compensation for yaw calculation from magnetometer and accelerometer
Assumptions

AFAIK, the sequence of the rotations are also important. It is not given in the question. I will assume the sequence to go from inertial frame to body frame is (Yaw, Pitch, Roll) = (Z, Y, X). So the sequence to go from body to inertial is the reverse.
Since these details are not present in the question, I assume the convention given in Link 1 and Link2.
Accelerometer frame is same as IMU reference

The transformation
The measured vectors are obtained in IMU (accelerometer) frame. To take a vector resolved in IMU frame, to inertial(?) frame, the transformation as given in the above reference is
$$
\begin{bmatrix}
v
\end{bmatrix}^I_{3\times1}
=
\begin{bmatrix}
    C\psi &amp; -S\psi &amp; 0\\
    S\psi &amp; C\psi &amp; 0\\
    0 &amp; 0 &amp; 1\\
\end{bmatrix}
\color{red}{
\begin{bmatrix}
    C\theta &amp; 0 &amp; S\theta\\
    0 &amp; 1 &amp; 0\\
    -S\theta &amp; 0 &amp; C\theta\\
\end{bmatrix}
\begin{bmatrix}
    1 &amp; 0 &amp; 0\\
    0 &amp; C\phi &amp; -S\phi\\
    0 &amp; S\phi &amp; C\phi\\
\end{bmatrix}
}
\begin{bmatrix}
v
\end{bmatrix}^{IMU}_{3\times1}
$$
The red matrices indicate what I assume is the transformation equation set shown in the question.
Assume the magnetometer data was available in the same frame of reference as the accelerometer. Let that reading be \$[x_M', y_M'z_M']^T\$.
$$
\begin{bmatrix}
x_{M2}\\
y_{M2}\\
z_{M2}
\end{bmatrix}
=
\begin{bmatrix}
C\theta &amp; S\phi S\theta &amp; C\phi S\theta\\
0 &amp; C\phi &amp; -S\phi\\
\dots &amp; \dots &amp; \dots
\end{bmatrix}
\begin{bmatrix}
x_{M}'\\
y_{M}'\\
z_{M}'
\end{bmatrix}
$$
Since the Y and Z axes are inverted for the magnetometer, the above equation changes to
$$
\begin{bmatrix}
x_{M2}\\
y_{M2}\\
z_{M2}
\end{bmatrix}
=
\begin{bmatrix}
C\theta &amp; S\phi S\theta &amp; C\phi S\theta\\
0 &amp; C\phi &amp; -S\phi\\
\dots &amp; \dots &amp; \dots
\end{bmatrix}
\begin{bmatrix}
x_{M}\\
\color{red}{-}y_{M}\\
\color{red}{-}z_{M}
\end{bmatrix}
$$
The above is significantly different from your equations.
Sanity check
You have mentioned in the comments that &quot;pitch on y axis&quot;. This means that a rotation about pitch should leave the Y component of a vector unchanged (if it was the last operation performed). Equation for yM in the question doesn't seem to satisfy that logic. Of course, This check is only correct assuming a certain sequence of rotations.
Note
I see that your equations seem almost correct if the sequence of rotations to go from inertial frame to body frame is (Yaw, Roll, Pitch). The negation on Y and Z components before applying the equations are still required. So, see if your output becomes correct if you insert
yM = -yM;
zM = -zM;

just before the transformation.