I use the model registration module to register bones and calculate the 3-d displacements.
My aim is to define the axis of rotation (finite helical axis).
I would like to calculate a mean rotation for all transforms (I can leave the translations for the minute) and then calculate the variability. (Ultimately to calculate the mean rotation axis and the variability around the mean axis)

An algorithm is provided for computing axis of rotation for a single pair of transforms and a perfect rotation around a single axis. Nobody volunteered to develop a more robust algorithm that could approximate the rotation axis from many imperfect measurements, so you would need to figure it out.

It may not be very hard, as there are not too many unknowns (you only have to determine position and orientation of the rotation axis and for each pair of transforms a rotation angle, so if you have 4 pair of transforms then in total you have 10 unknowns) and you may be able to get an approximate initial guess based on anatomical knowledge. You need to specify some error metric (such as residual position error of some anatomical landmark points) and then you can compute the uknowns using any non-linear optimizer.