Can we characterize the probability generating function as a linear operator? For a nonnegative integer-valued random variable $X$ with $\mathbb P(X=j)=p_j$, we define the probability generating function of (the distribution of) $X$ by $$P_X(s):=\mathbb E\left[s^X\right] = \sum_{j=0}^\infty p_js^j. $$ For a probability mass function $\{p_j\}$, define $\mathbb G(\{p_j\})=\mathbb E\left[s^X\right]$. Then the domain of $\mathbb G$ is $$D(\mathbb G)=\left\{\{p_j\}\in\mathbb R^{\mathbb N}: p_j\geqslant0,\ \sum_{j=0}^\infty p_j=1 \right\}. $$ Now, $D(\mathbb G)$ is not a vector space as it is not closed under addition. The most structure I have been able to find is that $D(\mathbb G)$ is that it is the $\sigma$-convex hull of the standard basis $\{e_n\}$ for $\ell^1(\mathbb R)$. 
For contrast, the characteristic function of (the distribution of) $X$ is well-defined as the Fourier-Stieltjes transform $$F_X\mapsto \int_{\mathbb R}e^{itx}\ \mathsf dF_X(x), $$ which is a bounded linear operator from $L^1(\mathbb R)$ to $L^\infty(\mathbb R)$.
The other question is how to characterize the image of $\mathbb G$ - all I've determined so far is that it is the $\sigma$-convex hull of $\{1,s,s^2,\ldots\}$ in the set of power series which converge absolutely on $[-1,1]$, which is contained in $C^\infty_c(-1,1)$. Meanwhile, e.g. Bochner's theorem tells us that characteristic functions of probability distributions correspond uniquely to normalized positive definite continuous functions on $\mathbb C$.
I realize the question is a bit open-ended so I have tagged it with (soft-question). I'm mostly looking for some intuition on what the probability generating function transform really is, in an analytic sense.
 A: Well one place we could start would be to note that we could define  on a larger set; namely any real sequence of coefficients for which the generating series is absolutely convergent or even Abel summable (namely the limit of the series exists as s→1 even if the series diverges when one substitutes s=1). If we are willing to consider quasiprobability distributions, then it seems like an extension of the probability generating function to some subspace of ℓ1 or $bs$ (see: en.wikipedia.org/wiki/Sequence_space#Other_sequence_spaces) would be fairly natural. 
I don't think intuitively that the sequences of coefficients such that probability generating series is absolutely convergent would correspond quite to all of ℓ1, but it would make the probability generating function better correspond to the Laplace transform (I think), since we're dropping the condition that the sequence must be in the convex cone.
This also means allowing negative $p_i$, which would remove some of the probabilistic meaning of the transform. However, the Laplace and Fourier transforms don't inherently have significance exclusive to probability theory either.
The Laplace and Fourier transforms are both defined as linear operators on larger spaces than just positive functions. If we want the probability generating function to be a linear operator, then it has to accept negative inputs -- otherwise if we limit it only to probability distributions we will be stuck to the convex cone, on which it is impossible to define a truly linear function (i.e. without restrictions on homogeneity). One reason why Laplace and Fourier transforms can be characterized as linear operators is because their definition is not restricted to a convex cone.
As it stands the probability generating function inherently has for its domain discrete random variables non-negative support, the way the Laplace transform has for its domain random variables with non-negative support. If we want a generalization similar to that of the Fourier transform (or double-sided Laplace transform) we would also need to consider discrete functions with possibly negative support. 
This perhaps would not be as difficult as it sounds, however, since this would correspond to the generalization of Taylor series to Laurent series, which is relatively unproblematic, although then we also need to consider the intricacies involved with domains of convergence being annuli instead of simply discs. Also we probably need to decide whether or not we want an analog of Fourier transform on actual functions or just generalized functions.
