Trouble with ICRS/TEME conversions

I’m trying to do some conversions between the ICRS and TEME reference frames. I need to get a direction cosine matrix (DCM) relating the two frames. The trouble I’m having is that the matrices I get are not orthogonal.

from astropy.coordinates import SkyCoord, ICRS, TEME
from astropy.time import Time
import numpy as np

def getMatrixICRSToTEME(time):
    return SkyCoord(x=np.array([1,0,0]),
                    y=np.array([0,1,0]),
                    z=np.array([0,0,1]),
                    representation_type='cartesian',
                    frame = ICRS()
                   ).transform_to(TEME(obstime=time)).cartesian.get_xyz().value

def getMatrixTEMEToICRS(time):
    return SkyCoord(x=np.array([1,0,0]),
                    y=np.array([0,1,0]),
                    z=np.array([0,0,1]),
                    representation_type='cartesian',
                    frame = TEME(obstime=time)
                   ).transform_to(ICRS()).cartesian.get_xyz().value

if __name__ == '__main__':
    time = Time("2024-05-28")
    icrs2teme = getMatrixICRSToTEME(time)
    print(icrs2teme @ icrs2teme.T)
    teme2icrs = getMatrixTEMEToICRS(time)
    print(teme2icrs @ teme2icrs.T)

This should result in identity matrices if the DCMs are orthogonal–and it’s close, but not to floating-point precision:

[[ 9.99999081e-01  5.35123344e-05  7.42733356e-05]
 [ 5.35123344e-05  1.00000058e+00 -5.15245529e-05]
 [ 7.42733356e-05 -5.15245529e-05  1.00000034e+00]]
[[ 9.99999074e-01 -5.39986969e-05 -7.41490657e-05]
 [-5.39986969e-05  1.00000059e+00  5.18876115e-05]
 [-7.41490657e-05  5.18876115e-05  1.00000034e+00]]

If I take the angle between the columns of the DCMs, I get about 0.003 degrees away from perpendicular:

x = icrs2teme[:,0]
y = icrs2teme[:,1]
z = icrs2teme[:,2]
print(np.rad2deg(np.arccos(np.dot(x,y))))
print(np.rad2deg(np.arccos(np.dot(x,z))))
print(np.rad2deg(np.arccos(np.dot(z,y))))

Returns:

89.99694037393135
89.9957605307974
90.00298212776633

Any insights would be appreciated, thank you!