I have the polynomial
$$f(x;a)=3ax^6+6x^5-9ax^4-4x^3+9ax^2+6x-3a$$
where the variable $a$ is a random variable from the uniform distribution in the range $[0,1)$. When I analyze this function using Mathematica and specific values for $a$, I find there are four complex roots, a negative real root, and a positive real root - I need to extract the positive real root to get an equation $x=x^+(a)$.
I have tried both of Kulkarni's methods${}^{1,2}$ in which he suggests splitting the sextic equation into a quadratic and a quartic or a cubic and another cubic, but numerics in Mathematica show this cannot be done.
At this point, I'm leaning toward using the Kampe de Feriet functions - these are referenced in many places but never given explicitly.
Could someone provide a hint to a solution, or a reference where I could find the application of the Kampe de Feriet function to sextic equations?
${}^1$http://elib.mi.sanu.ac.rs/files/journals/tm/21/tm1124.pdf ${}^2$http://euclid.trentu.ca/aejm/V3N1/Kulkarni.V3N1.pdf
Here's what the plot looks like near the region I'm interested in, for $a\in[0.1,0.9]$:

Actually, your sextic has a solvable Galois group so can be exactly expressed in terms of radicals. The solution turns out to be quite simple. Given
$$3a x ^6 + 6x^5 - 9a x^4 - 4x^3 + 9a x^2 + 6x - 3a = 0\tag1$$
Then the two real roots are,
$$x = u\pm\sqrt{u^2+1}\tag2$$
where,
$$u=\frac{1}{3a}\Big(z_1^{1/3}+\frac{1}{z_1^{1/3}}-1\Big)\tag3$$
and $z_1$ is any non-zero root of,
$$z^2+(9a^2+2)z+1=0\tag4$$
Your sextic's positive real root $x$ is the positive case $\pm$ of $(2)$.
Note: Be careful when testing this numerically, especially using $(3)$. What I find supremely annoying about my version of Mathematica is it gives a complex value to $z^{1/3}$ when $z$ is a negative real number.
P.S.: However, an easy way to verify this is to use resultants and the Mathematica syntax,
which then yields $(1)$, showing that your sextic in just a quadratic $(2)$ in $x$, with coefficients that is cubic $(3)$ in $u$.