Weak equivalence principle is not valid for stars - observations of mass-luminosity relation for visual binaries

Stars are special objects in the sense of the possible gravitational deviations: they have both barionic and non-barionic matter (like trapped and slowly advancing to the surface light) inside. While for barionic matter the weak equivalence principle was established by Kepler (planets orbiting around stars) it was never checked for stars.
The good example are binaries. There are many binary stars which are visible as double stars with resolved period and axis and ratio of inertial masses (through measurements of the velocities of stars). Many parameters of such stars are published in [1]
The usual formula applied to the stars from the third Kepler Law:

T²=4π²*a³/[G(m₁+m₂)]         (1)

Here T is the period of rotation of one star around the second one, a is semi-axis, m1 and m2 are masses of the stars (assuming gravitational mass is equal to inertial mass) and G is gravitational constant.
However, the light theoretically may have much higher gravitational pull compare to the inertial mass from E=mc*c relation (it is assumed that the inertial mass of light being emitted and reabsorbed inside star is still according to E=mc*c, as it was proved by Einstein himself) [2]. The presence of slow light may modify the gravitational pull, making it much stronger for the star which has more trapped light (and other non-barionic matter). While the exact amount of trapped light is difficult to calculate (not much is known about the light content of the interior of fully ionized plasma), it is obvious that this amount is correlated with luminosity of the star - the higher the luminosity, the higher the amount of trapped light and the higher the additional gravitational pull on the star (the higher the deviation between the gravitational and inertial mass).
In the derivation of the formula (1) the gravitational masses are always comes as a product [3]:

F=G*M1*M2/r²

Here M1 and M2 are gravitational masses. Assuming the added pull is proportional to luminosity which is proportional to mass (whether gravitational or inertial) [1], it is possible to assume:

F=G*K1*K2*m1*m2/r²

Here K1 and K2 are multiplicity coefficients, the value of K may be especially high to ultra-bright star. It is important that both coefficients for binaries are always a product.
The modified third Kepler Law:

T²=4π²*a³/[G*K1*K2*(m1+m2)]

Here m1 and m2 are inertial masses. When K1=K2=1, the third Kepler Law for barionic matter is obtained.
To determine the masses from the observation of binaries we need: T, a, and ratio of masses m1/m2=n. Since the ratio of masses is determined through the Doppler shift of spectra of stars, is a ratio of inertial masses. We have two equations for masses m1, m2:

G*K1*K2*(m1+m2)=4π²*a³/T²

m1/m2=n
Then:

m2=4π²*a³/[G*T²*K1*K2*(n+1)]
m1=4π²*a³*n/[G*T²*K1*K2*(n+1)]

Suppose we decided to determine the inertial masses from the visual binaries with two distinct masses m1>>m2. How it would influence the mass-luminosity correlation (like in [1])?
It is possible to show, that contrary to the case of valid third Kepler Law the slope of the dependence will be depended upon the ratio of masses!
 Lets  consider three cases:
1.Binary m1 and m1
2.Binary m2 and m2
3.Binary m1 and m2
In the first case the value of m1 is (because n=1)
m1=m1(old)/[K1*K1], here m1(old)=4π²*a³/[G*T²*2]
Here m1(old) is real inertial mass. K1 is large and the value of m1 is shifted strongly toward smaller mass compare to real inertial mass.
In the second case the value of m2 (n is equal to 1)
m2=m2(old)/[K2*K2]
If K2 is smaller than 1 (supposedly Sun has the value of K exactly one) the mass of smaller star will shifted strongly toward larger mass
In the third case the value of m1 is
m1=m1(old)/[K1*K2], m1(old)=m1(old)=4π²*a³*n/[G*T²*(n+1)]
Since both coefficients K1 and K2 are here, one is small and one is big, the shift compare to the real inertial mass is smaller
m2=m2(old)/[K1*K2]
This idea may be immediately checked. If the mass-luminosity curve is plotted using first only stars with close masses, it will be compressed because of K1*K1 and K2*K2 coefficients along the x-axis (the slope will be larger). If the same curve is plotted using the stars with different masses (preferably with large difference, but I used what we have in [1]) it will be much more spread. I manually chose approximately half of visual binaries (17 binaries or 34 stars) from Table 1A from [1] with close masses and obtained the relation between the luminosity and mass:
Absolute luminosity= -3.2119*ln(m)+5.1264
And for the rest of the binaries (21 binaries, 42 stars) - masses are different:
Absolute luminosity = -2.495*ln(m)+5.4042
The white dwarfs were excluded, like in [1].
The scattering in the second case is much larger (as expected, because the product of coefficients K1 and K2 is highly unpredictable), but they must give much smoother curve if the coefficients are the same - the shift is larger, but it is more predictable - obviously it is some smooth function of luminosity and can not jump from star to star).
Indeed as predicted the slope is larger beyond any error for the subset of close in masses stars compare to the far in masses stars.
Since the deviation from Newton law for the barionic matter of such large scale would be noticed long ago in our Solar system (some small deviation due to GR are not considered), the only possible explanation is that the weak equivalence principle does not hold for gravitation of stars (mixture of barionic and non-barionic matter), possibly due to mechanism outlined in [2].

References.
1.George E. McCluskey, Jr, Yoji Kondo "On the mass-luminosity relation" // Astrophysics and Space Science 17 (1972), p.p.134-149
http://articles.adsabs.harvard.edu//full/1972Ap%26SS..17..134M/0000137.000.html
2.D.S.Tipikin, publication on blog
https://tipikin.blogspot.com/2019/10/stars-are-full-of-trapped-light-may.html
3.http://www.astro.caltech.edu/~george/ay20/Ay20-Lec4x.pdf