It might interest you to know that you can represent a 4D rotation by a pair of quaternions.

This is related to a matrix form by isoclinic decomposition.

The representation is two unit quaternions qL, qR and the transformation of 4D vector v is

qL * v * qR

(quaternion multiplication)

I did some 4D visualization software using this as the rotation representation, I still have it hanging around but got tired of it

See

https://en.wikipedia.org/wiki/Rotations_in_4-dimensional_Euclidean_space