Let us note $^A\vec\gamma$ a vector of magnitude $|\vec\gamma|=\gamma$ defined in the frame $A$. Let $A$ and $B$ be two distinct frames, $|^A\vec\gamma| = |^B\vec\gamma|$. Moreover, it may be useful to recall the definition of an inertial frame: it is a frame whose acceleration is nil.
Finally, let's recall the transport theorem, as per "Analytic mechanics of space systems" by Schaub and Junkins (this is an exact quote):
Let $N$ and $B$ be two frames with a relative angular velocity vector of $\vec\omega_{B/N}$, and let $\vec r$ be a generic vector; then the derivative of $\vec r$ in the $N$ frame can be related to the derivative of $\vec r$ in the $B$ frame as: $$ \frac{^Nd}{dt} \vec r= \frac{^Bd}{dt} \vec r + \vec\omega_{B/N} \times \vec r$$
A subsequent quote of great interest is the following:
... we will find that vectors are typically differentiated with respect to an inertial frame called $N$.
As explained in section 2.4 of "GPS" by G. Xu and Y. Xu, 2016:
The motion of satellites follows Newtonian mechanics, [and these] are only valid and expressed in an inertial coordinate system.
In the case of astrodynamics, the words "frame" and "coordinate system" are almost always interchangeable.
We can now answer your questions:
Is it mandatory to convert to ECEF? I would be happy to get a detailed explanation on how the orientation affects the acceleration. The propagation of the motion of the spacecraft due to external gravitional forces must be computed in an inertial frame. If the motion of that spacecraft is mostly due to Earth's gravity, then the celestial motion must be computed in the ECI frame. However, it is mandatory to convert the satellite's position and velocity to the ECEF frame to compute how the spherical harmonics of Earth affect the satellite's motion. For an in-depth definition of the ECEF frame, I highly recommend Chapter 2 of "GPS" by G. Xu and Y. Xu, 2016 (I believe this chapter is freely available on the Springer website). For a further explanation of why this is important, note that if one frame is inertial but not the other, then the magnitude of a given vector at a given time may be different, but the components of said vector won't be. For example, imagine two frames $A$ and $B$ which, at time $t_k$ are oriented exactly such that the conversion from $A$ to $B$ corresponds to a rotation of $+\pi/2$ about the Z axis. Further, let $^A\vec\alpha = [1,~0,~0]^T$. Then, in the $B$ frame, we have $^B\vec\alpha = [0,~-1,~0]^T$. Therefore, if an acceleration in the $A$ frame leads to a change of the $x$ component of $-1$ at a subsequent time of $t_n$, $^A\vec\alpha_{t_n} = [0,~0,~0]^T$. But it would be wrong to apply that change from the frame $A$ directly to $^B\vec\alpha$ without transforming that vector into the $A$ frame first. In fact, if you would, you'd find that your (incorrect) updated $^B\vec\alpha_{t_n} = [-1,~-1,~0]^T$. That's evidently wrong since $|^A\vec\alpha_{t_n}| \neq |^B\vec\alpha_{t_n}|$.
How to convert the acceleration to a different coordinate system? I can convert the coordinates but not the change in velocity. You don't convert the acceleration to a different frame actually. You convert the position and velocity to the ECEF frame, compute the acceleration due to the harmonics in ECEF, and then rebuild the updated position and velocity in the ECEF frame. Finally, you convert your updated ECEF state back into ECI. As a proof, this is how it's done in GMAT 2016a. (I'm linking to Github because the source code is significantly easier to browse there.)
I hope this helps.
