In[]:=
Clear["Global`*"];(****CreateDataforTesting*******)nn2=2;nbv=7;rr=10;data0=Table[{x1,x2,x3,x4,x5,x6,x7,1.2x1-3.4x2+10x3+2.8x4+1.6x5-1.9x6+2.8x7+RandomReal[rr]},{x1,RandomReal[rr,nn2]},{x2,RandomReal[rr,nn2]},{x3,RandomReal[rr,nn2]},{x4,RandomReal[rr,nn2]},{x5,RandomReal[rr,nn2]},{x6,RandomReal[rr,nn2]},{x7,RandomReal[rr,nn2]}];data0=Flatten[data0,nbv-1];data=Drop[data0,None,{nbv+1}];resp=Drop[data0,None,{1,nbv}];resp=Flatten[resp];(******Fitamodellm1*******)
In[]:=
lm1=LinearModelFit[{data,resp}]s1=(lm1["FitResiduals"])^2;
(*******Fitthesquarederrorstoamodelwiththesamepredictors\lm2****)
In[]:=
lm2=LinearModelFit[{data,s1}];
(******ComputetheBreusch–Paganteststatistic*****)
bp=(Variance[s1](Length[data]-1)-Total[lm2["FitResiduals"]^2])/2/(Total[s1]/Length[data])^2(******OR*******)
bp2=With[{sqResids=lm1["FitResiduals"]^2},(Variance[sqResids](Length[data]-1)-Total[lm2["FitResiduals"]^2])/2/(Total[sqResids]/Length[data])^2]
(*********Computethep-value****)
dr=Length[lm1["BestFitParameters"]]-1;1-CDF[ChiSquareDistribution[dr],bp]