I found a nice problem recently, but could not come up with a solution:
Find all functions $f:\mathbb{Z}_{\geq 0}\to\mathbb{Z}_{\geq 0}$ such that for all $0< \lambda < 1$, if $X\sim G(\lambda)$, then there exists $0<\mu < 1$ with $f(X) \sim G(\mu)$.
Thanks to the answers below, I was able to understand how to solve it!
Spoiler: the solutions are all the functions $f_d:x\mapsto \lfloor x/d\rfloor$ for a fixed $d\geq 1$.
Let $f$ be a function such that for each $\lambda\in (0,1)$, if $X\sim\mathcal G(\lambda)$, then there exists $\mu\in (0,1)$ such that $f(X)\sim\mathcal G(\mu)$.
First of all, $f$ should satisfy $f\left(\mathbb Z_{\geqslant 0}\right)=\mathbb Z_{\geqslant 0}$ because $\mathbb P(f(X)=k)=\mathbb P(f^{-1}(\{k\}))$ is positive for each $k$ and $X$ satisfying a geometric distribution.
For $j,k\in\mathbb Z_{\geqslant 0}$, define $a_{k,j}=\begin{cases}1&\mbox{if }f(j)=k\\ 0&\mbox{otherwise} \end{cases}$, $a_{k,-1}=0$ and $I_k=f^{-1}(\{k\})$.
We observe that if $X\sim\mathcal G(\lambda)$, then $f(X)\sim\mathcal G(\mu(\lambda))$ where $$ \mu(\lambda)=\sum_{j\in I_0}(1-\lambda)^j\lambda. $$ Indeed, $\mu(\lambda)=\mathbb P(f(X)=0)=\mathbb P(X\in I_0)$.
It will be more convenient to write $\mu(1-\lambda)$ as a power series, namely, $$ \mu(1-\lambda)=\sum_{j\in I_0}\lambda^j-\sum_{j\in I_0}\lambda^{j+1} =\sum_{j=0}^{\infty}\left(a_{0,j}-a_{0,j-1}\right)\lambda^j. $$
By definition of geometric distribution with parameter $\mu(1-\lambda)$, we derive that for each $k\geqslant 0$ and each $\lambda\in (0,1)$,
$$\tag{*} \left(1-\sum_{j=0}^{\infty}\left(a_{0,j}-a_{0,j-1}\right)\lambda^j\right)^k \sum_{j=0}^{\infty}\left(a_{0,j}-a_{0,j-1}\right)\lambda^j =\sum_{j=0}^{\infty}\left(a_{k,j}-a_{k,j-1}\right)\lambda^j. $$ Note that all the involved series are convergent since $0\leqslant a_{k,j}\leqslant 1$ and $0<\lambda<1$.
Let us look first at the constant coefficient, which can be obtained by letting $\lambda\to 0$. We get that for each $k\geqslant 0$, $$ \left(1-a_{0,0}\right)^ka_{0,0}=a_{k,0} $$ If $a_{0,0}=0$, we would get that $a_{k,0}=0$ for each $k\geqslant 0$ but since the sets $I_k$ are disjoint and their union is $\mathbb Z_{\geqslant 0}$, we should have $\sum_{k\geqslant 0}a_{k,0}=1$ hence we get $a_{0,0}=1$ hence $f(0)=0$.
Let us now update $(*)$. Separating the index $j=0$ in the two sums of the left hand side, we derive that for each $k\geqslant 1$ and each $\lambda\in (0,1)$,
$$\tag{**} \left(-\sum_{j=1}^{\infty}\left(a_{0,j}-a_{0,j-1}\right)\lambda^j\right)^k \left(1+\sum_{j=1}^{\infty}\left(a_{0,j}-a_{0,j-1}\right)\lambda^j\right) =\sum_{j=0}^{\infty}\left(a_{k,j}-a_{k,j-1}\right)\lambda^j. $$
Notice that in the left hand side of $(**)$, the lowest degree in $\lambda$ is $k$ hence we should have for $0\leqslant j\leqslant k-1$ that $a_{k,j}-a_{k,j-1}=0$, and since $a_{j,-1}=0$, $a_{k,0}=a_{k,1}=\dots=a_{k,k-1}=0$ or in other words, $j<k$ implies that $f(j)\neq k$ hence $f(j)\leqslant j$.
Looking now at the coefficient in front of $\lambda^k$, we get that $(1-a_{0,1})^k=a_{k,k}$. If $f(1)=1$, then $f(k)=k$ for all $k$.
Suppose that $f(1)\neq 1$. By 6, we have $f(1)=0$. Let $d$ be the smallest natural number such that $f(d)\neq 0$, that is, $a_{0,d}=0$ and $a_{0,j }=1$ for $j<d$. Looking now at the coefficient in front of $\lambda^j$, $0\leqslant j\leqslant j_0k$, we find that $f(dk)=k$.
Define $b_{k,j}=a_{k,j}-a_{k,j-1}$ and $c_j$ by $c_0=0$ and $c_j=-b_{0,j}$ for $j \geqslant 1$. Since $-\sum_{j=1}^{\infty}\left(a_{0,j}-a_{0,j-1}\right)\lambda^j=\sum_{j=0}^\infty c_j\lambda^j$, for $k=1$ in $(**)$, recognizing a Cauchy product between power series, we get that \begin{align} b_{1,j}&=\sum_{i=0}^j b_{0,i}c_{j-i}\\ &=-\sum_{i=0}^{j-1} b_{0,i}b_{0,j-i}\\ &=-b_{0,j}-\sum_{i=1}^{j-1} b_{0,i}b_{0,j-i}\\ &=-b_{0,j}-\sum_{i=1}^{j-1} b_{0,i}\mathbf{1}_{i\geqslant d}b_{0,j-i}\mathbf{1}_{j-i\geqslant d} \end{align} where we used in the last equality that $b_{0,i}=0$ if $1\leqslant i<d$. We thus got $$\tag{1} b_{1,j}=-b_{0,j}-\mathbf{1}_{j\geqslant 2d}\sum_{i=d}^{j-d} b_{0,i}b_{0,j-i} $$
We will now prove by induction on $\ell\geqslant 1$ that $b_{0,d+\ell}=0$.
Since $$ b_{1,2d+1}+b_{0,2d+1}=\mathbf{1}_{f(2d+1)\in\{0,1\}}-\mathbf{1}_{f(2d)\in\{0,1\}}\in \{-1,0,1\}, $$ this forces $b_{0,d+1}$ to be equal to $0$.
Going back to (**), we can drastically simplify the left hand side by $$ \left(-\sum_{j=1}^{d}\left(a_{0,j}-a_{0,j-1}\right)\lambda^j\right)^k \left(1+\sum_{j=1}^{d}\left(a_{0,j}-a_{0,j-1}\right)\lambda^j\right) =\sum_{j=0}^{\infty}\left(a_{k,j}-a_{k,j-1}\right)\lambda^j. $$ and since $a_{0,j}-a_{j-1}=0$ for $1\leqslant j\leqslant d-1$, we get $$ \lambda^{dk} \left(1-\lambda^d\right) =\sum_{j=0}^{\infty}\left(a_{k,j}-a_{k,j-1}\right)\lambda^j $$ hence identifying the coefficient in $dk$ and $d(k+1)$ gives that $f(k)=\lfloor k/d\rfloor$, where $\lfloor x\rfloor$ is the unique integer for which $\lfloor x\rfloor\leqslant x<\lfloor x\rfloor+1$.
It remains to check that any function of the form $f(k)=\lfloor k/d\rfloor $ for $d\geqslant 1$ does the job. First, $I_0=\{0,\dots,d-1\}$ hence $\mu(\lambda)=1-(1-\lambda)^{d}$ and if $\ell$ is a non-negative integer, \begin{align} \mathbb P(f(X)=\ell)&=\mathbb P(d\ell\leqslant X\leqslant d\ell+d-1)\\ &=\sum_{j=d\ell}^{d\ell+d-1}(1-\lambda)^j\lambda&\\ &=(1-\lambda)^{d\ell}\sum_{j=0}^{ d-1}(1-\lambda)^j\lambda&\\ &=(1-\lambda)^{d\ell}\left(1-(1-\lambda)^{d}\right) \\ &=(1-\mu(\lambda))^{ \ell}\mu(\lambda)&. \end{align}