Assuming the ESD structure activates in the same way (a supply-shorting mechanism seems to imply we're talking clamp diodes here), you could reproduce the ESD waveform, scaled appropriately. 4kV HBM implies 2.67A peak and a time constant of 0.15us. So, set Rs and Cmax accordingly.
For C higher than the scaled value, you can't guarantee the same pulse duration, and so more charge will be absorbed, which might activate other mechanisms (e.g. electromigration). Energy or power may be above their respective curves as well, but you don't know those curves, how they scale with time that is, so the next safest plan is to honor the DC current limit.
So, from available data that would be something like a limit of 2.67A for some 10s of ns, dropping along some slope or curve to the 10mA DC figure. Almost certainly, this curve/slope can be shallower than the discharge tail of the HBM event, but taking that as a minimum limit, it drops to 10mA after about 800ns, so that seems a safe assumption.
Example: so, whereas HBM is 100pF + 1500R at 4kV, you could presumably do, say, 1nF at 400V and 150R. If supply was 15V, maybe even as much as 27nF with 5.6R. Mind, the pin's internal resistance will start to dominate in this regime, dissipating more out of the total energy (but, 27nF @ 15V is also much less energy than 100pF @ 4kV, so that's not necessarily a problem), and also affecting the time constant -- that is to say, we can't assume the pin is an ideal ground down here, whereas if it's dropping, you know, 30V or something, that's negligible out of 4kV.
The real problem, then, is trying to infer scaling laws from limited data, and extrapolating across two decades of scale!
In any case, the given problem is trivially avoided with external clamp diodes (which will handle much more current at lower voltage drop), so that's not a big deal. Supply filtering may be relevant as well (assuming the shorting event happens elsewhere; and if it does happen locally, it probably means something local has failed worse to begin with?).
Tim