Let $$ f(x) = \sum_{k=1}^\infty\frac{(-1)^{k+1}x^k}{(k-1)!\zeta(2k)}$$
Are there lower bounds, upper bounds or (unlikely) simpler closed form for $f(x)$?
The bounds for Bernoulli numbers I found give the trivial bounds.
Wolfram Alpha times out.
mpmath's nsum returns small negative values.
The Wikipedia page on the Riesz function which Barry Cipra cited has a link to a paper by Wilf showing that there are infinitely many real $x$ where the function vanishes, though rather sparse: the number of such $x \geq X$ is asymptotically proportional to $\log X$. The function certainly does not approach zero as $x \rightarrow \infty$; even the estimate $f(x) = o(x^{1/4+\epsilon})$ would be equivalent to the Riemann hypothesis, as Riesz noted in his original 1916 paper in Acta Math.
Numerically, the first few zero crossings beyond $1.15$ are at about $19326.551$, $22521.798$, $51868.607$ if I did this right. Here's a plot for $x \in [5\cdot 10^3, 10^5]$, comparing $f(x)$ with the main term $$ {\rm Re}\left( \frac {\Gamma\left(1 - \frac{\rho_1}{2}\right)}{\zeta'(\rho_1)} x^{\rho_1/2} \right) $$ in the asymptotic expansion of $f(x)$ for large $x$.
Here $\rho_1 = \frac12 + 14.1347\ldots i$ is the first complex zero of the Riemann $\zeta$ function. The factor $x^{\rho_1/2}$ oscillates with amplitude $x^{1/4}$. The coefficient is small (absolute value about $7.775 \cdot 10^{-5}$) because of the complex Gamma factor; further complex zeros provide additional terms, but the $\rho_2$ term is already smaller by a factor of almost $300$, and further terms should be smaller yet, assuming no $\zeta$ zeros get very close to each other (let alone collide or worse). The discrepancy between the two plots is accounted for almost entirely by the first trivial zero of $\zeta$ at $-2$, which gives a term $1/\left(2\zeta'(-2)x\right)$ that is asymptotically negligible but makes $f(x)$ negative for $1.15 < x < 5000$. It takes a long time for the growth of $x^{1/4}$ to overcome the small Gamma factor: $\left|\phantom.f(x)\right|$ first exceeds $1$ around $x = 10^{19}$.
The asymptotic expansion comes from the contour-integral formula $$ f(x) = \frac1{2\pi i} \int_{c-i\infty}^{c+i\infty} \frac{\Gamma(s+1)}{\zeta(-2s)} x^{-s} ds $$ with $-1 < c < -1/2$. NB the denominator really is $\zeta(-2s)$, not $\zeta(+2s)$; we encounter zeros of $\zeta$ by moving the contour to the right. The Wikipedia page gives this formula but without the factor of $(2\pi i)^{-1}$ $-$ maybe somebody reading this can fix that error.
To compute $f(x)$ numerically, we can't quite use the formula $$ f(x) = \sum_{n=1}^\infty \mu(n) \frac{x}{n^2} e^{-x/n^2} $$ because it doesn't converge fast enough to locate the large zeros. Instead, by applying the same trick to the sum over (say) $k \geq 4$ instead of $k \geq 1$ we get an expansion such as $$ f(x) = \frac{6x}{\pi^2} - \frac{90x^2}{\pi^4} + \frac{945x^3}{2\pi^6} + \sum_{n=1}^\infty \mu(n) \frac{x}{n^2} \left( \exp\left(\frac{-x}{n^2}\right) - 1 + \frac{x}{n^2} - \frac{x^2}{2n^4} \right) $$ which was used to compute $f(x)$ to enough accuracy to draw the plot shown above. For large $x$ this technique still requires $x^{1/2 + \epsilon}$ terms. It may be barely feasible to compute $f$ this way for $x \approx 10^{19}$ large enough to see $\left|\phantom.f(x)\right|>1$; but the asymptotic expansion using the first handful of zeta zeros should be much better for this purpose.
This answer is inspired by a now deleted answer by Abhinav Kumar:
A perhaps simpler form is $$\sum_{n=1}^{\infty} \mu(n) \frac{x}{n^2} e^{- x/n^2}$$ where $\mu$ is the Möbius function.
Proof: Expand $1/\zeta(2k) = \sum_{n=1}^{\infty} \mu(n) n^{-2k}$ and rearrange (everything is absolutely convergent) to get $$\sum_{k=1}^{\infty} \frac{(-1)^{k+1} x^k}{(k-1)!} \sum_{n=1}^{\infty} \frac{\mu(n)}{n^{2k}} = \sum_{n=1}^{\infty} \mu(n) \sum_{k=1}^{\infty} \frac{(-1)^{k-1}}{(k-1)!} \left( \frac{x}{n^2} \right)^{k} = \sum_{n=1}^{\infty} \mu(n) \left( \frac{x}{n^2} \right) \exp \left( \frac{- x}{n^2} \right). \quad \square$$
For $x$ small, we get an easy bound of $\pi^2/6 x$ by replacing $\mu(n)$ by $1$ and $\exp(-x/n^2)$ by $1$. It should be possible to push this down to $(6/\pi^2+o(1)_{x \to 0}) x$. Cut the sum at some fixed $N$. Then we have $$\sum_{n=1}^{N} \mu(n) \frac{x}{n^2} \exp(-x/n^2) = \sum_{n=1}^N \left( x \frac{\mu(n)}{n^2} (1 - O(x/n^2) \right) = \sum_{n=1}^{N} \left( x\frac{\mu(n)}{n^2} + O(x^2/n^4) \right) = x (6/\pi^2 - O(1/N)) + O(x^2)$$ and $$\sum_{n=N+1}^{\infty} \mu(n) \frac{x}{n^2} \exp(-x/n^2) = O\left(x \sum_{n=N+1}^{\infty} \frac{1}{n^2} \right) = O(x/N)$$
For large $x$, the problem seems much harder. If we look at the part of the sum where $0.9 \leq x/n^2 \leq 1.1$, we are basically getting $e^{-1} \sum_{n=\sqrt{x}/\sqrt{1.1}}^{\sqrt{x}/\sqrt{0.9}} \mu(n)$. My understanding is that sums like $\sum_{n=ay}^{by} \mu(n)$ for fixed $a$ and $b$ are expected to behave like $c \sqrt{y}$, so this seems to be contributing $c \sqrt[4]{x}$, contrary to the numerical observation that it is going to $0$. I can't imagine any easy method getting a bound which approaches $0$, in light of this observation.
UPDATE. Barry Cipra points out that this is the Riesz function. It is known to have infinitely many real zeroes, contrary to the numeric evidence, and the estimate $f(z) = O(z^{1/4+\epsilon})$ is equivalent to the Riemann Hypothesis (which makes sense, because I pointed out that it is related to bounds like $\sum_{n=ay}^{by} \mu(n) = O(\sqrt{y})$, which are related to RH.)
