How to Decompose and Compose the local transform matrix The Setup:
I have a simple setup with two cube objects with a parent child relationship like so:

The Red cube is a child object of the green cube.
The Red cube has the following transforms:

The green cube has the following transforms:

The Problem
The thing I am confused about is, when I try to decompose the local_matrix of the red cube and compose it back to its original transform matrix. I get different results. Below is a script I used to decompose and  compose the local matrix. Any clues why they are different? And what is the correct way to re-compose the local matrix given the translation, rotation and scale matrices.
Script and Result
import bpy
import mathutils
import math

myObject = bpy.data.objects['Cube.002']

print(&quot;\n location mat&quot;)
loc = myObject.matrix_local.to_translation()
print(loc)
mat_loc = mathutils.Matrix.Translation(loc)
print(mat_loc)

print(&quot;\n rotation mat&quot;)
rot = myObject.matrix_local.to_quaternion()
print(rot)
mat_rot = rot.to_matrix().to_4x4()
print(mat_rot)

print(&quot;\n scale mat&quot;)
scl = myObject.matrix_local.to_scale()
print(scl)
mat_scl = mathutils.Matrix.Identity(4)
mat_scl_x = mathutils.Matrix.Scale(scl.x, 4, (1,0,0))
mat_scl_y = mathutils.Matrix.Scale(scl.y, 4, (0,1,0))
mat_scl_z = mathutils.Matrix.Scale(scl.z, 4, (0,0,1))
mat_scl = mat_scl_x @ mat_scl_y @ mat_scl_z
print(mat_scl)

print(&quot;\n myLocalMat mat&quot;)
myLocalMat =  mat_loc @ mat_rot @ mat_scl
print(myLocalMat)

print(&quot;\n local mat&quot;)
print(myObject.matrix_local)

 location mat
&lt;Vector (0.0000, 0.0000, 3.0000)&gt;
&lt;Matrix 4x4 (1.0000, 0.0000, 0.0000, 0.0000)
            (0.0000, 1.0000, 0.0000, 0.0000)
            (0.0000, 0.0000, 1.0000, 3.0000)
            (0.0000, 0.0000, 0.0000, 1.0000)&gt;

 rotation mat
&lt;Quaternion (w=0.9280, x=0.0000, y=0.0000, z=0.3726)&gt;
&lt;Matrix 4x4 (0.7224, -0.6915, 0.0000, 0.0000)
            (0.6915,  0.7224, 0.0000, 0.0000)
            (0.0000,  0.0000, 1.0000, 0.0000)
            (0.0000,  0.0000, 0.0000, 1.0000)&gt;

 scale mat
&lt;Vector (1.5811, 1.5811, 1.0000)&gt;
&lt;Matrix 4x4 (1.5811, 0.0000, 0.0000, 0.0000)
            (0.0000, 1.5811, 0.0000, 0.0000)
            (0.0000, 0.0000, 1.0000, 0.0000)
            (0.0000, 0.0000, 0.0000, 1.0000)&gt;

 myLocalMat mat
&lt;Matrix 4x4 (1.1422, -1.0934, 0.0000, 0.0000)
            (1.0934,  1.1422, 0.0000, 0.0000)
            (0.0000,  0.0000, 1.0000, 3.0000)
            (0.0000,  0.0000, 0.0000, 1.0000)&gt;

 local mat
&lt;Matrix 4x4 (1.4142, -1.4142, 0.0000, 0.0000)
            (0.7071,  0.7071, 0.0000, 0.0000)
            (0.0000,  0.0000, 1.0000, 3.0000)
            (0.0000,  0.0000, 0.0000, 1.0000)&gt;

Best Regards.
QC
Using the 3x3 matrix rotation part.
Does a child object inherit the matrix from the parent?
Can decompose a matrix via Matrix.decompose()
loc, rot, scale = M.decompose()

the rotation is given as a quaternion.  There is no guarantee that that the same rotation matrix will be returned when converting that quaternion back to a matrix.  https://en.wikipedia.org/wiki/Rotation_formalisms_in_three_dimensions#Rotation_matrix_%E2%86%94_quaternion
Test script, with both &quot;standard&quot; quaternion decompose and using the 3x3 rotation matrix part.
import bpy
from mathutils import Matrix

ob = bpy.context.object
# basis (or local) matrix
Mb = ob.matrix_local

print(&quot;Input&quot;)
print(Mb)


t, q, s = Mb.decompose()

T = Matrix.Translation(t)

R = q.to_matrix().to_4x4()

S = Matrix.Diagonal(s.to_4d())

M = T @ R @ S
print(&quot;Quaternion decompose&quot;)
print(M)

print((M.to_3x3() - Mb.to_3x3()).determinant())

R = Mb.to_3x3().normalized().to_4x4()

T = Matrix.Translation(
        Mb.to_translation()
        )

S = Matrix.Diagonal(
        Mb.to_scale().to_4d()
        )
        
M = T @ R @ S

print(&quot;3 x 3 decompose&quot;)
print(M)

print((M.to_3x3() - Mb.to_3x3()).determinant())

Sample output.
Input
&lt;Matrix 4x4 (1.8334, -0.5595,  0.0113, 0.6324)
            (0.5749,  0.1069, -0.9128, 2.6868)
            (6.2012,  1.4492,  0.3085, 0.8749)
            (0.0000,  0.0000,  0.0000, 1.0000)&gt;
Quaternion decompose
&lt;Matrix 4x4 ( 4.5809, -1.0592, -0.1913, 0.6324)
            (-1.7931, -0.0093, -0.9261, 2.6868)
            ( 4.2364,  1.1414, -0.1851, 0.8749)
            ( 0.0000,  0.0000,  0.0000, 1.0000)&gt;
0.6159369945526123
3 x 3 decompose
&lt;Matrix 4x4 (1.8334, -0.5595,  0.0113, 0.6324)
            (0.5749,  0.1069, -0.9128, 2.6868)
            (6.2012,  1.4492,  0.3085, 0.8749)
            (0.0000,  0.0000,  0.0000, 1.0000)&gt;
0.0