How to calculate bootstrap SE of the between-group difference in covariance (between two vars)

[shem shen]

Dear list,

My example code is as follows:

webuse auto, clear
forvalues i=1/50 {
    qui gen rhsbrw`i' = .
    qui rhsbsample, strata(foreign) cluster(make) weight(rhsbrw`i')
    qui replace rhsbrw`i'=rhsbrw`i'*weight 
}
qui svyset make [pw=weight], strata(foreign) vce(bootstrap) bsrweight(rhsbrw*) mse
corr headroom gear_ratio, cov

(obs=74)

             | headroom gear_r~o
-------------+------------------
    headroom |  .715707
  gear_ratio | -.145857  .208198

svy: "how to calculate bootstrap SE of the diff in covariance between headroom and gear_ratio"?

In this example, suppose I have two separate samples (marked by the foreign variable), and both samples have clustered structures (clustervar is make). Further, I have a sampling weight variable weight. How should I use the bootstrap methods to derive the standard error of the difference in the covariance between headroom and gear_ratio? I provided some code above using Philippe Van Kerm's rhsbsample command. However, I am stuck at the last step. One may say that I can simply treat covariance as the mean of the product [X-E(X)][Y-E(Y)], and do the following:

sum head
local headmn=r(mean)
sum gear
local gearmn=r(mean)
gen covv=(head-`headmn')*(gear-`gearmn')
svy: reg covv i.foreign

However, it should be noted that the point estimate of the covariance is not simply the mean of [X-E(X)][Y-E(Y)]. The denominator is weighted N - 1, not weighted N. This makes me believe that I cannot simply use the above syntax for my purpose. See below:

. sum covv

    Variable |        Obs        Mean    Std. Dev.       Min        Max
-------------+---------------------------------------------------------
        covv |         74    -.143886    .3080931   -1.30679   .4317239

. corr head gear ,cov
(obs=74)

             | headroom gear_r~o
-------------+------------------
    headroom |  .715707
  gear_ratio | -.145857  .208198

The mean of covv is actually slightly different from what we can get from corr, cov, because of the difference in the denominator.
Thank you!

Reference:
Philippe Van Kerm, 2013. "Repeated half-sample bootstrap resampling," United Kingdom Stata Users' Group Meetings 2013 10, Stata Users Group.