Solving quintic equations of the form $x^5-x+A=0$

The answer to your third question is yes! The method uses Bring radicals, whose explicit form in terms of generalized hypergeometric functions can be found using the Lagrange inversion theorem. (In fact since any quintic can be reduced to this form, in principle this method can be used to solve any quintic.) I can answer your second and fourth questions partially by developing this method. But I'm afraid I will only be able to obtain the first solution with coefficients $(A, 0, 0, 0)$.

The idea at its core is quite simple. Basically we rewrite the equation as $x^5 - x = - A$, treat the left hand side as a function $f(x) = x^5 - x$, then try to answer the question "what is $f^{-1}(-A)$." This is then done by expressing $f^{-1}$ as a power series. The Lagrange inversion theorem gives this inverse as $$ x = \sum_{k=0}^\infty \binom{5 k}{k} \frac{A^{4k+1}}{4k+1}\ . $$ Unfortunately, this series doesn't converge for all values of $A$. In fact the radius of convergence is $4/(5\times 5^{1/4})\approx 0.535$, so evaluating the series directly would only give us the solution for one integer $A = 0$. This is where the generalized hypergeometric function comes in. We can analytically continue this series to define a function of $A$. The function whose power series (at zero) is $$ \sum_{n=0}^{\infty}\prod_{k=0}^{n} \frac{(k+a_1)\cdots(k+a_p)}{(k+b_1)\dots(k+b_q)(k+1)} z $$ is denoted as $_p F_q(a_1,\dots, a_p;b_1,\dots,b_q;z)$. To convert our function $f^{-1}(A)$ into the standard form, we need to compute the ratio between consecutive terms, which is $$ \begin{align} & \quad \frac{(5k +5)!A^{4k+5}}{(k+1)!(4k+4)!(4k+5)}\cdot\frac{k!(4k)!(4k+1)}{(5k)!A^{4k+1}}\\ & = \frac{(5k+5)(5k+4)(5k+3)(5k+2)(5k+1)(4k+1)A^4}{(k+1)(4k+4)(4k+3)(4k+2)(4k+1)(4k+5)} \\ & = \frac{5(5k+4)(5k+3)(5k+2)(5k+1)}{4(4k+5)(4k+3)(4k+2)(k+1)}A^4 \\ & = \frac{(k+1/5)(k+2/5)(k+3/5)(k+4/5)}{(k+1/2)(k+3/4)(k+5/4)(k+1)}\left(5\left(\frac{5A}{4}\right)^4\right)\ . \end{align} $$ Now since the numerator has four factors and the denominator has three factors besides $(k+1)$, this is $_4F_3$ (times an extra factor of $A$, since the starting term in our series is $A$, not $1$). The parameters are the numbers added to $k$ in each factor, and the argument is $(5^5/4^4)A^4 = (3125/256)A^4$. This gives you the first solution $A \;_3F_4(\frac{1}{5}, \frac{2}{5}, \frac{3}{5}, \frac{4}{5}; \frac{1}{2}, \frac{3}{4}, \frac{5}{4}; \frac{3125}{256}A^4)$.

In order to obtain the other roots, in principle we can find use this root to factor the polynomial and try to solve the resulting quartic. However that involves too much computation and doesn't seem like the neatest way of obtaining the results you got here.

I have a full answer to your question now, and it turns out to be rather long and indirect. There are multiple ways of arriving at these solutions to the quintic $x^5 - x + t= 0$, but in this answer I will only pursue one route. It involves solving a differential equation associated to the original quintic, called the differential resolvent. This idea originally came from Cockle (1860) and Harley (1862), but I also drew from this question on the Mathematica Stack Exhange. I should note that this different method due to M. L. Glasser possibly takes fewer pages to explain and is in some ways more natural, but it relies on the Lagrange reversion theorem (not to be confused with the Lagrange inversion theorem), which is also not particularly straightforward to apply.

Let me first describe the patterns that you observed a bit more systematically, then I will dive into the actual solution process, whence the patterns will emerge at various stages.

The patterns

The four generalized hypergeometric functions are related. (I'll call them GHFs from now on.) You probably already noticed that their arguments are the same. The parameters of $_4F_3\left(\frac{-1}{20},\frac{3}{20},\frac{7}{20},\frac{11}{20};\frac{1}{4},\frac{1}{2},\frac{3}{4};\frac{3125 t^4}{256}\right)$ are in fact roots of a quartic (for the first group) and a cubic (for the second group) that will come up in the solution process. For notational convenience, let's call the first group of parameters $a_1, a_2, a_3, a_4$ and the second group $b_1, b_2, b_3$. The other three GHFs are obtained by replacing each $a_i$ and $b_i$ by $1+a_i - b_j$ and $1+b_i-b_j$ for $i = 1, 2, 3$, expect for $b_j$ itselft which becomes $2-b_j$. For example, when $j = 1$, we have $b_j = \frac{1}{4}$, so we have $$ \begin{align} \left(\frac{-1}{20},\frac{3}{20},\frac{7}{20},\frac{11}{20}\right) & \mapsto \left(1+\frac{-1}{20}-\frac{1}{4},1+\frac{3}{20}-\frac{1}{4},1+\frac{7}{20}-\frac{1}{4},1+\frac{11}{20}-\frac{1}{4}\right) \\ & = \left(\frac{7}{10},\frac{9}{10},\frac{11}{10},\frac{13}{10}\right) \end{align} $$ and $$ \begin{align} \left(\frac{1}{4},\frac{1}{2},\frac{3}{4}\right) & \mapsto \left(2 - \frac{1}{4}, 1+\frac{1}{2}-\frac{1}{4}, 1 + \frac{3}{4} - \frac{1}{4}\right)\\ & = \left(\frac{7}{4}, \frac{5}{4}, \frac{3}{2}\right)\ . \end{align} $$ Thus we get the second GHF in your expression. (Note that the order of GHF parameters don't matter.)

The coefficients $B_1, B_2, B_3, B_4$ will turn out to be the first four terms in a series that expresses $x$ in terms of $t$. This is why you see consecutive powers of $t$ in them. Curiously, even though they will be the last thing we look for, they can actually be evaluated after the first step.

The derivation

Finding the differential resolvent

This approach starts essentially the same as the previous one — we take the equation $x^5 - x + t = 0$ and treat $x$ as a function of $t$. Instead of directly finding the power series of $x(t)$, however, Cockle and Harley's method constructs a differential equation for $x(t)$, called the differential resolvent, which is solved using GHFs. The differential resolvent is a fourth order ordinary differential equation of the following form, $$ \mu_0(t) x''''(t) + \mu_1(t) x'''(t) + \mu_2(t) x''(t) + \mu_3(t) x'(t) + \mu_4(t) x(t) + \mu_5 = 0\ , $$ where $\mu_i(t)$ are polynomials in $t$. You may be wondering why the equation is fourth order, which would only have four linearly independent solutions, when we are looking for five solutions. The reason is, since the $x^4$ term in $x^5 - x + t = 0$ is zero, the roots should always add up to zero, so there should be at most four linearly independent solutions.

Constructing the differential resolvent is a three-step process. First we implicitly differentiate the original quintic repeatedly and iteratively solve for the derivatives of $x$ up to the fourth order. For example, differentiating once we get $5 x^4 x' - x' + 1 = 0$, so $x' = 1/(1-5x^4)$. Differentiating twice gives us $20 x^3 (x')^2 + 5 x^4 x'' - x'' = 0$, and since we already know $x'$ in terms of $x$, we plug that in and find that $x'' = 20x^3/(1-5x^4)^3$. Continuing this way, we can express all derivatives of $x$ in terms of $x$.

Next we plug these expressions into the differential resolvent and get a polynomial equation in $x$, $t$, and $\mu_i$. This equation has up to $x^29$ in it, but fortunately we can use the original quintic to reduce that. Basically we replace every $x^5$ by $x-t$, or equivalently find the remainder of the polynomial division by $x^5-x+t$. In the end the highest remaining power of $x$ in the expression should be less than $5$.

The third step is to set the coefficient of each power of $x$ to zero. This give us five equations which we can solve for $\mu_i$ in terms of $t$. These equations are all linear in $\mu_i$, so they are definitely manageable. Note that there are six unknowns and five equations, so the extra degree of freedom allows us to choose $\mu_i$'s to all be polynomials in $t$. In the end, (one possible form of) the differential resolvent for this quintic is $$ (3125 t^4 - 256) x'''' + 31250 t^3 x''' + 73125 x'' + 31875 x' - 1155 x = 0\ . $$

We expect a solution to the quintic to also solve the differential resolvent, because we have essentially used the quintic to construct a linear combination of derivatives of $x$ that is congruent to zero modulo $x^5 - x + t$.

Solving the differential resolvent

Using the knowledge (or one might call it premonition :P) that one of the solutions to this equation is of the form $_4F_3(a_1, a_2, a_3, a_4;b_1, b_2, b_3;3125 t^4/256)$, we seek to cast this as a generalized hypergeometric-type differential equation, $$ \begin{multline} z \frac{d}{dz} \left(z \frac{d}{dz} + b_1 - 1\right) \left(z \frac{d}{dz} + b_2 - 1\right) \left(z \frac{d}{dz} + b_3 - 1\right) x(z)\\ = z \left(z \frac{d}{dz} + a_1\right) \left(z \frac{d}{dz} + a_2\right) \left(z \frac{d}{dz} + a_3\right) \left(z \frac{d}{dz} + a_4\right) x(z)\ . \end{multline} $$ It is known that the four linearly independent solutions are $_4 F_3(a_1, a_2, a_3, a_4; b_1, b_2, b_3; z)$ along with three other GHFs whose parameters are derived using the method described earlier, where the GHF derived from $b_i$ is multiplied by a prefactor $z^{1-b_i}$.

Thus we would like to show that the differential resolvent can be obtained by substituting $3125 t^4/256$ for $z$ and expanding all the derivatives. So we do that, then match up the coefficients of each derivative of $x(t)$ and each power of $t$, and we get four equations for $a_p$ and three for $b_q$, which are up to quartic in $a_p$ and up to cubic in $b_q$. Luckily they can all be rewritten as linear equations of elementary symmetric polynomials in $a_p$ and $b_q$. (This shouldn't come as a surprise, after all the GH equation is invariant under permutation of $a_p$ and $b_q$.) We can easily solve for these symmetric polynomials, which are exactly the coefficients of the quartic (resp. cubic) that has $a_p$ (resp. $b_q$) as its solutions. The resulting quartic and cubic turn out to be pretty easy to solve, and we get $$ \begin{gather} a_1 = - \frac{1}{20}\ ,\quad a_2 = \frac{3}{20}\ ,\quad a_3 = \frac{7}{20}\ ,\quad a_4 = \frac{11}{20}\ ,\\ b_1 = \frac{1}{4}\ ,\quad a_2 = \frac{1}{2}\ ,\quad b_3 = \frac{3}{4}\ . \end{gather} $$ From this we can obtain the four linearly independent solutions to the differential resolvent, $$ \begin{align} F_1 & = \;_4F_3\left(\frac{-1}{20},\frac{3}{20},\frac{7}{20},\frac{11}{20};\frac{1}{4},\frac{1}{2},\frac{3}{4};\frac{3125 t^4}{256}\right)\\ F_2 & = \frac{125}{64} 5^{3/4} t^3 \;_4F_3\left(\frac{7}{10},\frac{9}{10},\frac{11}{10},\frac{13}{10};\frac{5}{4},\frac{3}{2},\frac{7}{4};\frac{3125 t^4}{256}\right)\\ F_3 & = \frac{25\sqrt{5}}{16} t^2 \;_4F_3\left(\frac{9}{20},\frac{13}{20},\frac{17}{20},\frac{21}{20};\frac{3}{4},\frac{5}{4},\frac{3}{2};\frac{3125 t^4}{256}\right)\\ F_4 & = \frac{5}{4} 5^{1/4} t \;_4F_3\left(\frac{1}{5},\frac{2}{5},\frac{3}{5},\frac{4}{5};\frac{1}{2},\frac{3}{4},\frac{5}{4};\frac{3125 t^4}{256}\right) \end{align} $$

We reasoned that a solution to the original quintic should solve the differential resolvent, and thus take the form $x = c_1 F_1 + c_2 F_2 + c_3 F_3 + c_4 F_4$. However, for an arbitrary set of coefficients $c_i$, this general solution will satisfy the differential resolvent but not necessarily the original quintic. In principle, in order to figure out the right coefficients, we need to substitute the general solution into the original quintic and solve the resulting equations, which seems like a harder problem than what we set out to tackle! Fortunately, there is a shortcut. We could check that the series expansion of the solution at $t = 0$ satisfies the equation. Recall that each GHF is of the form $\sum_{k=0}^\infty \alpha_k t^{4k}$, where $\alpha_0 = 1$, so the first four terms of the series are exactly $B_4, B_1, B_2$, and $B_3$. This tells us that $(B_4, B_1, B_2, B_3) = (x(0), x'(0) t, x''(0) t^2/2, x'''(0) t^3/6)$. Since we already know how to write the derivatives of $x$ in terms of $x(t)$, we can simply plug the five solutions to $x^5-x=0$ into the expressions we used in step one of the derivation of the differential resolvent to evaluate these coefficients. Finally, we find that $$ (B_4, B_1, B_2, B_3) = (0, t, 0, 0)\ , $$ or $$ (B_4, B_1, B_2, B_3) = \left(\omega, -\frac{t}{4}, - \frac{5 \omega^3 t^2}{32}, -\frac{5\omega^2 t^3}{32}\right)\ , $$ where $\omega$ are the four quartic roots of unity.

This concludes the procedure, and we arrive at the solution you found using Mathematica. Hopefully this also demonstrates the origin of the patterns you observed.