library(lmtest)
library(car)
getwd()
setwd("C:/Users/花寻觅/Documents/xuexuexue")
a=read.csv("C:/Users/花寻觅/Documents/xuexuexue/四川异方差.csv")
h<-lm(Y~X,data=a)
summary(h)
第一种检验方法：残插图检验
#plot(h,which=1)which还可以等于2,3,4
e2<- residuals(h)^2
x1<- a$X
x2<-a$T
par(mfrow=c(1,2))
plot(x1,e2,main="残差图")
plot(x2,e2,main="残差图")
最后观察，但是数据足够大才好做

第二种检验方法：BP检验#Breusch-Pagan
代码：
bptest(h，studentize=FALSE)去掉学生化修正异方差的作用，如果p-value小于0.05可以认为存在异方差
#BPTEST原理summary(lm(resid(h)^2)~X,data=a)
Goldfeld-Quandt test：GQ检验的思想是对数据回归得到残差序列,然后把残差作为被解释变量,原方程各解释变量作为解释变量做回归,得到bp的统计量

看p值， <0.05 拒绝原假设 得到存在 异方差

             >0.05不能拒绝原假设 即回归不存在异方差代码：gqtest(h)
第三种检验方法：WHITE
bptest(h,~fitted(h)+I(fitted(h)^2))#WHITE检验
如果是两个解释变量，则有：white<- lm(e2~X*T+I(T^2)+I(X^2)+X+T,data=a)#省略也可以
summary(white)
用R2*样本容量与查表（0.05显著水平下的已知自由度），如果大于则拒绝同方差假设
解决措施：
####加权回归，权数为：weights=1/ (a$ X)^2，
#(h1<-glm(Y~X,data=a,weights=1/ (a$ X)^2))？这里我也疑惑为什么和summary(h)运行出来一样
#标准误法
coeftest(h)
coeftest(h,vcov=hccm)