Dear Rick,

Thank you for this very nice introduction into Passing-Bablok regression and a possible implementation in SAS!
Like the late Wolfgang Bablok, I work at Roche Diagnostics, Germany. Therefore, the use and discussion of Passing-Bablok regression has a strong tradition in our company. Some of the insights gained since have recently been published (see below).

But let me first make some remarks on the bootstrap estimator for the CIs of the slope:
The effect of the duplicate values on the bootstrap estimator of the variance or CI is certainly not as high as as one might fear following your argumentation. Using the Poisson bootstrap, it is easy to see why. In the Poisson bootstrap one does not draw exactly samples of size N, but only with an average sample size N with poissonian distribution. Then, the probability to draw each of the original samples is again Poissonian with mean 1. The ratio of samples dropped is 1/e=37%, as you stated. You can also calculate the average number of pairs between equal replicates which is $ 1/e \sum_{i=2}^\infty 1/(2*(i-2)!)=0.5 $ per sample. So from the N(N-1)/2 pairs, only N/2 are between replicates in a bootstrap sample and therefore their effect is asymptotically negligible. (Btw., also the least squares estimator can be formulated in terms of pairwise differences, so that a similar influence of repeated values is to be expected for the LS estimator.) Obviously using other types of bootstrap where exactly N samples are drawn will not change this conclusion significantly.
The asymptotical confidence intervals given by Passing-Bablok may be too narrow, at least for strongly heteroscedastic data. Recently, I proposed to use a more general expression, which uses the fact that the Passing-Bablok estimator can be derived from finding the slope value where a certain Kendall correlation coefficient vanishes:
https://www.degruyter.com/document/doi/10.1515/ijb-2019-0157/html?lang=de
Kendall himself derived a more general expression for the variance of his coefficient, which can be used to derive CI's for the PB slope, cf. equation 4 in this preprint I wrote together with Prof. Jakob Raymaekers:
https://arxiv.org/pdf/2202.08060.pdf
At least to leading order in N, bootstrapping seems to recover this estimator.

After their first paper, Wolfgang Bablok and Heinrich Passing tried to extend their method to situations when the expected slope is far from 1. They dubbed this situation "method transformation" in contrast to "method comparison", see https://edoc.hu-berlin.de/bitstream/handle/18452/11842/cclm.1988.26.11.783.pdf?sequence=1
Recently, Terry Therneau, in the vignette of his "deming" R-package (https://CRAN.R-project.org/package=deming), could show, that in one of these methods, the equivariant PB regression, the slope estimator is simply the median of the absolute values of all pairwise slopes, which is significantly more easy to calculate than the usual PB estimator, which may best be considered as a first order approximation to this more general estimator. (The sign of the slope follows from the sign of Kendalls correlation coefficients of the measurement pairs).
In my paper cited above, this estimator is derived from a Kendall correlation and it is shown that it is a robust version of a method known as either "reduced major axis regression" or "geometric mean regression". A fast quasilinear implementation is available in the "robslopes" 'R-package by Jakob Raymaekers (cf. the preprint and https://CRAN.R-project.org/package=robslopes) and will be included in a new version of our "mcr" package.
The traditional PB regression estimate is probably best considered as an approximation to the equivariant estimator. Therefore, I recommend to give the equivariant PB regression a try and hope, that it soon will replace the older variant and will find a place in SAS, too!

Kind regards
Florian Dufey