OK, this is more interesting than I thought (and I'd welcome comments/criticism from others). Salvatore's question is prompted by a difference in results between gretl's Breusch-Pagan test for heteroskedasticity (via the command modtest --breusch-pagan) and the test as implemented in Eviews. First point: on the face of it there's a difference between the test actually proposed by Breusch and Pagan (Econometrica, 1979), which I'll call Test1, and that performed in Eviews, which I'll call Test2. (Eviews is by no means alone in using Test2, it's also the one described in the Wikipedia entry.) Both tests involve auxiliary regressions estimated via OLS, with the same set of independent variables, but * In Test1 the dependent variable is the squared residual from the original regression divided by the ML estimator of the error variance for that regression, and the LM test statistic is half of the explained sum of squares from the auxiliary regression. (This is very clearly set out on the first page of the B-P paper.) * In Test2 the dependent variable is the (unscaled) squared residual from the original regression and the LM test statistic is sample size times R^2 from the auxiliary regression. At first I thought this difference might be merely apparent (that is, the results should be numerically identical if done right). But I couldn't prove that and (unless I've messed up) simulation shows that they're not at all numerically identical. (Maybe they're asymptotically equivalent?) In that case, which one is "right" or "better"? Monte Carlo to the rescue. I show my test script below. It can be run in two modes, depending on the value of "H0_true" near the top of the script. * If H0_true = 1 the simulated error is homoskedastic and the rejection rates for the two tests give a measure of their sizes. I'm using the 5 percent significance level, so ideally the empirical rejection rate should be close to 0.05. * If H0_true = 0 the simulated error is by construction heteroskedastic and the rejection rates give a measure of power. On running the script in both modes several times, here's what I see: * Size: both tests seem about right. Sometimes one is closer to the nominal size of 0.05, sometimes the other. * Power: Test1 (the original) seems to be substantially better. Here are rejection rates from a few runs: H0_true = 1 Test1 0.049, Test2 0.053 Test1 0.043, Test2 0.046 Test1 0.048, Test2 0.047 H0_true = 0 Test1 0.548, Test2 0.311 Test1 0.546, Test2 0.316 Test1 0.698, Test2 0.520 If anyone sees an error in my simulation, please let me know! <hansl> set verbose off # sample size for regressions nulldata 100 # uncomment below for reproducible results # set seed 76543 # error term is homoskedastic? H0_true = 1 # independent variable series x = normal() # 5 percent critical value for chi-square(1) scalar X1crit = critical(X, 1, 0.05) # number of replications K = 5000 # recorders matrix reject = zeros(1,2) matrix LMdiff = zeros(K,1) loop i=1..K --quiet e = normal() if H0_true # homoskedastic y = 10 + x + e else # heteroskedastic y = 10 + x + 0.2*x*e endif # estimate original regression ols y 0 x --quiet u2 = $uhat^2 # Test1: as Breusch and Pagan (1979) vml = $ess / $T g = u2/vml ols g 0 x --quiet LM1 = (sst(g) - $ess)/2 reject[1] += LM1 > X1crit # Test2: as Wikipedia, Eviews ols u2 0 x --quiet LM2 = $trsq reject[2] += LM2 > X1crit LMdiff[i] = LM1 - LM2 endloop reject /= K printf "Rejection rates at 5 percent level:\n" printf " Test1 %.3f, Test2 %.3f\n", reject[1], reject[2] colnames(LMdiff, "LMdiff") summary --matrix=LMdiff --simple </hansl> Allin

On Sun, 14 Jun 2020, Sven Schreiber wrote: Certainly relevant; thanks, Sven. I hadn't twigged that "Test2" is equivalent to Koenker's robust B-P version. I re-ran my test script with uniform errors (should probably try some other cases) and found: * Under H0 the original B-P test is "under-sized": rejects at much less than 5 percent frequency using alpha = 0.05. The size of the robust version is roughly right. * Under my H1, error=0.2*x*uniform(), the original B-P test still rejects with much higher frequency than the Koenker variant. Allin

On Mon, 15 Jun 2020, Sven Schreiber wrote: I tried another run of my Monte Carlo, with t(12) errors. Original B-P was somewhat oversized (around 0.075 or 0.08 for nominal 0.05) but reasonably powerful. The Koenker robust variant was correctly sized at 0.05 (as claimed) but had about 1/2 to 5/8 the power against my H1 (with error multiplied by 0.2*x). It seems the subsequent literature has kind of skated over Koenker's 1981 caveat, that the statistic is "not entirely satisfactory" since "the power of the resulting test may be quite poor except under idealized Gaussian conditions". I would add, _even_ under Gaussian conditions. Allin