Is there an analytic way of solving univariate polynomial equations in general? (All polynomials here are assumed to be univariate)
I know that elementary solutions in terms of the coefficients may not exist for polynomials of order 5 or higher, but Wikipedia says that quintic and sextic equations could be done in analytic form using special functions such as the Jacobi theta.
My question is- do such analytic (not necessarily elementary) solutions of even higher polynomial equations also exist?
 A: More a long comment than an answer. The solution to the problem of finding an analytic expression for the roots of the general polynomial equation of order $n$
$$
a_n z^n + a_{n-1} z^{n-1} + \ldots + a_1 z + a_0=0 \label{1}\tag{1}
$$
was, according to Giuseppe Belardinelli ([1] pp. 3-4), triggered by the discovery of hypergeometric functions of $n$ variables. Precisely, by using such classes of functions, two analytic formulas, one by Hjalmar Mellin and one by Giuseppe Belardinelli himself, were found in 1921: a brief description of their approaches is given below.
Mellin (see [2] and also [1], §20, where the transformation of \eqref{1} in \eqref{1a} is explicitly described) starts by considering an equivalent form of equation \eqref{1}, precisely the following one
$$
z^n+b_1 z^{n_1}+\ldots+b_mz^{n_m} - 1=0,\label{1a}\tag{1a}
$$
where $n_i<n$ and $b_j\neq 0$ for each $j=1,\ldots,m<n-1$. This form has the property that if one of its solution solutions is known, say $z_p=z_p(b_1,\ldots,b_m)$, it is possible to find all the remaining ones by the following formula
$$
z_i(b_1,\ldots,b_m)=\zeta_iz_p(\zeta_i^{n_1}b_1,\ldots,\zeta_i^{n_p}b_m)\quad i=1,\ldots, n\label{2}\tag{2}
$$
where $\zeta_i$, $i=1,\ldots, n$ are the solutions of the cyclotomic equation
$$
\zeta^n-1=0
$$
Then Mellin succeeds in finding a particular hypergeometric function $z_p$, which he calls solution principale, such that $z_i(b_1,\ldots,b_m)$ is a root of \eqref{1a}, and thus by \eqref{2} he obtains a closed form analytic expression for the roots of \eqref{1}.
The solution of Belardinelli is similar in that it relies again on hypergeometric functions, but it is conceptually different: following the footsteps of Alfredo Capelli, he seeks the $n$ roots $z_i$, $i=1,\ldots, n$ of \eqref{1} in the form of convergent power series: and he succeeds in this task by defining a multiple Lagrange series of $n$-th order Pochammer hypergeometric functions that gives exactly the roots of \eqref{1} as a function of the coefficients $a_0,\ldots,a_n$ for each finite $n$ (see [1] for the details).
Additional note
The only English language reference I was able to find on the same topic is the short (translated from the Russian) paper [3] of E. N. Mikhalkin, which builds his result on the one of Mellin. He gives a very simple form for $z_p(b)$ (theorem 1, formula 4 at p. 302 of [3]): however, his formula is under all aspects a simplification of the one given by Mellin's in reference [2], and he doesn't explain nor describes how Mellin gets its result, nor provides any historical survey. 
Bibliography
[1] Giuseppe Belardinelli (1960), Fonctions hypergéométriques de plusieurs variables et résolution analytique des équations algébriques générales (French), Mémorial des sciences mathématiques, no. 145 (1960), 80 p., MR121518, Zbl 0097.05901.
[2] Hjalmar Mellin (1921), Résolution de l’équation algébrique générale à l’aide de la fonction gamma (French), Comptes Rendus Hebdomadaires des Séances de l’Académie des Sciences, Paris, 172, pp. 658-661, JFM 48.1238.02.
[3] Evgenii Nikolaevich Mikhalkin, On solving general algebraic equations by integrals of elementary functions, Siberian Mathematical Journal 47, No. 2, 301-306 (2006). MR2227983, ZBL1115.33001.
