I'll leave alone the OpenPGP spec, and consider the problem of identifying among $k$ public keys $(n_i,e_i)$ the RSA public key $(n_j,e_j)$ used to encrypt $m$ messages per RSA with proper encryption padding. The only way is by examining the $m$ cryptograms $c_\ell$, which essentially are indistinguishable from uniformly random in $[0,n_j)$ (I assume the $c_\ell$ are computed with no countermeasure to foil identification attempts, except in last section).
An identification algorithm can use two facts:
- If $n_i\le\max c_\ell$ , that rules out $i=j$.
- The mean $x=\frac1m\sum c_\ell$ , which the attacker should compute, has Bates distribution, with mean $(n_j-1)/2$, standard deviation $(n_j-1)/\sqrt{12m}$. That distribution quickly becomes close to normal when $m$ increases, and the 68-95-99.7% rule quickly becomes applicable.
I think the optimum attacker's strategy in practice is to pick $j$ as $i$ with $n_i/2$ closest to $x$ among the $n_i$ surviving test 1. It's possible to estimate the confidence in a guess by computing $t_i=\lvert n_i/2-x\rvert\sqrt{12m}/n_i$ for the values of $i$ with $n_i/2$ the closest to $x$: the higher $t_i$, the less likely $i$ is (ignoring information obtained from other $n_i$): $t_i>1$ with 32% probability, $t_i>2$ with 5% probability, $t_i>3$ with 0.3% probability.
To go deeper in the analysis, we need an idea of the distribution of the $n_i$.
By making the say 256 high-order bits of the $n_i$ equal at generation time (fixed to some haphazard constant), which is easy and does not harm RSA security or interoperability, we could make identification hopeless.
For a common generation method with the prime factors $p_i$ and $q_i$ of $n_i$ independent and roughly uniform in $[2^{(b-1)/2},2^{b/2}]$ where $b$ is the bit size of $n_i$ (e.g. $b=3072$), I think an attacker would need $m$ in the order of a few times $k^2$ for a good success rate, and significantly more to conclude with good success rate and confidence for most users. Keys with the highest $n_j$ tend to be easier to guess correctly, because rule 1 prunes best, and there tends to be less other keys with nearby $n_i$. The second effect also makes keys with the lowest $n_j$ easier to guess than those in the center of the pack: those with the high byte near 0xB9 (for $b$ multiple of 8) tend to be the most resistant to guessing.
All this is amenable to simulation.
There's a possible countermeasure to foul identification attempts even if we don't want to change pre-existing $n_i$, : users of the system decide a common parameter $t$ corresponding to how many times it's tolerable to slow down the modular exponentiation part of encryption (e.g. $t=128$, which still won't make encryption more costly than decryption) and when computing a $c_l$ keep the smallest of $\left\lceil t\,n_i/2^b\right\rceil$ computations with different random padding.