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 in [0xB0,0xC1] (for $b$ multiple of 8) tend to be the most resistant to guessing, while giving no strong sign they have been chosen for this purpose.
All this is amenable to simulation and detailed analysis. I called for help. In the desert so far.
Update: this comment by Squeamish Ossifrage points a relevant paper: Jean Paul Degabriele, Victoria Fehr, Marc Fischlin, Tommaso Gagliardoni, Felix Günther, Giorgia Azzurra Marson, Arno Mittelbach and Kenneth G. Paterson: Unpicking PLAID - A Cryptographic Analysis of an ISO-standards-track Authentication Protocol, in proceedings of SSR 2014. Start with section 3.1 Tracing Cards via Shill-Key Ciphertexts.
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.