Subtract quaternions in python?
Rotation Difference.
Believe what you are looking for is rotation_difference ie the minimum rotation required to rotate from one quaternion to another
Python console code
Define a couple of angles using Eulers,
&gt;&gt;&gt; r1 = Euler((radians(30), radians(60), radians(45)))
&gt;&gt;&gt; r2 = Euler((radians(0), radians(60), radians(45)))

As quaternions
&gt;&gt;&gt; q1 = r1.to_quaternion()
&gt;&gt;&gt; q2 = r2.to_quaternion()

The difference
&gt;&gt;&gt; qd = q1.rotation_difference(q2)


&gt;&gt;&gt; list(map(degrees, qd.to_euler()))
[-30.000000834826057, 2.481181138556882e-06, -3.316504514323493e-06]

Or other way
&gt;&gt;&gt; qd = q2.rotation_difference(q1)
&gt;&gt;&gt; list(map(degrees, qd.to_euler()))
[30.00000424992064, -3.4388197927083468e-06, 1.2883001108333913e-06]

And using this result to get back to q1 from q2
&gt;&gt;&gt; qr = q2 @ qd
&gt;&gt;&gt; qr
Quaternion((0.8223630785942078, 0.022260043770074844, 0.5319756269454956, 0.20056211948394775))

&gt;&gt;&gt; q1
Quaternion((0.8223631381988525, 0.022260019555687904, 0.5319756865501404, 0.20056211948394775))

difference (within precision) is  pretty much Identity Quat (1, 0, 0, 0)
&gt;&gt;&gt; qr.rotation_difference(q1)
Quaternion((1.0, -1.4901161193847656e-08, 3.4458935260772705e-08, -2.8870999813079834e-08))