r_stuio回归函数的程序包是哪个?

Python014

r_stuio回归函数的程序包是哪个?,第1张

1、线性模型~回归分析:【包】:stats 【函数】:lm(formula, data, ...)逐步回归:step(lm(formula, data, ...))回归诊断:influence.measure(lm(formula, data, ...))多重共线性:kappa(XX,exact=T), eigen(XX)自相关检验:一阶:dwtest(y~x) 多阶:bgtest(y~x,order=2,type=”Chisq”)【备注】:1)stats包里的lm()可做多元线形模型,anova.mlm()比较多个多元线形模型,manova()做多元方差分析(MANOVA)。2)sn包的msn.mle()和 and mst.mle()可拟合多元偏正态和偏t分布模型。3)pls包提供偏最小二乘回归(PLSR)和主成分回归;4)ppls包可做惩罚偏最小二乘回归;5)dr包提供降维回归方法,如:片逆回归法(Sliced Inverse Regression)、片平均方差估计(sliced average variance estimation)。6)plsgenomics包做基于偏最小二乘回归的基因组分析。7)relaimpo包可评估回归参数的相对重要性。2、logistic回归:【包】:stats 【函数】:glm(formula, family=gaussian,data, ...)注:familybinomial(link = "logit") gaussian(link = "identity") Gamma(link = "inverse") inverse.gaussian(link = "1/mu^2") poisson(link = "log") quasi(link = "identity", variance = "constant") quasibinomial(link = "logit") quasipoisson(link = "log")

3、无监督分类~决策树:【包】:rpart 【函数】:rpart(formula,data, method="class",control=ct,parms=list(prior=c(p,1-p),split="information"))

rpart.plot(fit,branch=1,branch.type=2,type=1,extra=102,shadow.col=”gray”,box.col=”green”,

split.cex=1.2,main=”Kyphosis决策树”) #提供了复杂度损失修剪的修剪方法

printcp(fit):告诉分裂到哪一层,CP,nsplit,rel,error,交叉验证的估计误差(xerror),标准误差(xstd)

prune(fit,cp=fit$cptable[which.min(fit$cptable[,"xerror"]),"CP"]):剪枝函数

【备注】:1)CRAN的 MachineLearning任务列表有对树方法的细节描述。

2)分类树也常常是重要的多元方法,rpart包正是这样的包,

3)rpart.permutation包还可以做rpart()模型的置换(permutation)检验。

4)TWIX包的树可以外部剪枝。

5)hier.part包分割多元数据集的方差。

6)mvpart包可做多元回归树,

7)party包实现了递归分割(recursive partitioning),

8)rrp包实现了随机递归分割。

9)caret包可做分类和回归训练,进而caretLSF包实现了并行处理。

10)kknn包的k-近 邻法可用于回归,也可用于分类。

4、支持向量机:

【包】:e1071,kernlab

【函数】:svm(x_train,y_train,type="C-classification",cost=10,kernel="radial",probability=TRUE,scale=FALSE)

svp=ksvm(x,y,type="C-svc",kernel="rbf",kpar=list(sigma=1),C=1)

5、无监督分类~聚类分析:

【包】:stats

【函数】:系统聚类:hclust(d,method=”complete”,members=NULL)

快速聚类:kmeans(x,centers,iter.max=10,nstart=1,algorithm=“Hartigan-Wong”)

距离函数:dist(x,method=”euclidean”,diag=FALSE,upper=FALSE,p=2)

【备注】:1)CRAN的Cluster任务列表全面的综述了R实现的聚类方法。

2)stats里提供等级聚类hclust()和k-均值聚类kmeans()。

3)cluster包里有大量的聚类和可视化技 术,

4)clv包里则有一些聚类确认程序,

5)e1071包的classAgreement()可计算Rand index比较两种分类结果。

6)Trimmed k-means聚类分析可由trimcluster包实现,

7)聚类融合方法(Cluster Ensembles)由clue包实现,

8)clusterSim包能帮助选择最佳的聚类,

9)hybridHclust包提供一些混合聚类方法。

10)energy包里有基于E统计量的距离测度函数edist()和等级聚类方法hclust.energy()。

11)LLAhclust包提供基于似然(likelihood linkage)方法的聚类,也有评定聚类结果的指标。

12)fpc包里有基于Mahalanobis距离的聚类。

13)clustvarsel包有多种基于模型的聚类。

14)模糊聚类(fuzzy clustering)可在cluster包和hopach包里实现。

15)Kohonen包提供用于高维谱(spectra)或模式(pattern)的有监督和无监督的SOM算法。

16)clusterGeneration包帮助模拟聚类。

17)CRAN的Environmetrics任务列表里也有相关的聚类算法的综述。

18)mclust包实现了基于模型的聚类,

19)MFDA包实现了功能数据的基于模型的聚类。

地理加权回归(GWR)在R里面怎么实现?

121 人关注0 条评论

写回答

查看全部 5 个回答

写回答

叶山Shan Ye

GIS/地质/人文地理/可持续发展

A2A 谢邀,

我和我认识的一些人,刚开始用R做空间分析的时候,也遇到过这个问题。R这种开源的东西,优点是各种包很丰富,缺点是有些包的说明写得很乱,地理加权回归(GWR)的R包其实功能很强大,但大部分说明都不大靠谱。

GWR在R里面可以用好几个不同的包来实现,其中步骤最简单的是spgwr。思路就两步:建立窗口、用窗口扫全局。这其实就是GWR本质上的两步。比如我要在全美国范围内统计某两个(或多个)变量之间的回归关系,我可以做一个全局回归(global regression),但因为这些变量在空间分布上或许会有异质性(heterogeneity),表现在统计结果上就是空间不稳定性(nonstationarity),因此只看全局的统计,可能看不出什么结果来。举个不完全恰当但是很容易领会精神的例子,你比如说,我要分析亚洲范围内,经济发展程度与牛肉销量之间的关系,经济越发达的地方,人们就越吃得起牛肉。可是等我统计到印度的时候,坏了,印度大部分人不吃牛肉,这不是经济状况导致的,这一下就影响了全局统计的参考价值,那怎么办呢?我们可以建立一个窗口(正规说法是带宽窗口,bandwidth window),每次只统计窗口范围内的经济与牛肉销量的关系,然后用这个窗口去扫过全局的范围。等统计到印度的时候,印度内部的各地和印度自己比,吃牛肉的人的比例就不会突然减少,这样就能减少这种空间不稳定性对全局统计的影响。

所以,第一步就是要建立这样一个『窗口』。当然了,首先要安装包,我们要用到的R包有:

library(spgwr)

library(rgdal)

library(sf)

library(spData)

library(sp)

library(lattice)

library(ggplot2)

library(ggthemes)

其中,spgwr是做GWR的包,rgdal是用来读取矢量要素的,sf,sp和spData都是用来处理矢量数据的,别的基本都是画图用。

以下默认你会R和GWR的基本操作。并且,以下只展现方法,不要纠结我的数据和结果,我随便找的数据,这个数据本身没有什么意义,所以做出的统计看起来很『壮观』。

我们先导入数据。这里我用的是美国本土48州各个县(county,也有翻译成郡的)的人口普查数据和农业数据,来源是ESRI Online数据库。为啥用这个数据呢?因为...我电脑里面就存了这么个可以用来做GWR的数据...

我们用rgdal读取数据,然后把它画出来看看

require(rgdal)

usa_agri <- readOGR(dsn = "~/Documents/Spatial", layer = "usa_counties")

plot(usa_agri)

会得到这个东西:

readOGR里面,dsn后面加储存shp的路径(加到文件夹为止),layer后面写shp的文件名(不加.shp)。不喜欢rgdal的同学可以不用,用maptools或者spData等别的处理shp的R包代替。不过如果用maptools,要注意处理一下参考系。

我们看一下这个shp里面的列联表都有什么:

可见,shp里面有3108个县的数据,数据有61种。然后再看data下面有什么:

总之就是各种人口普查的数据,后面截不完图,还有经济、房地产和农业之类的数据。那我们就随便选两个来当变量。我就随便挑了,因变量选AVESIZE12,即2012年各个县农场的平均占地面积。自变量选POP_SQMI,也就是人口密度(每平方英里的人口)。

现在正式建立窗口,调用的是spgwr里面的gwr.sel函数:

bw <- gwr.sel( AVE_SIZE12 ~ POP_SQMI, data = usa_agri, gweight = gwr.Gauss,

verbose = FALSE, method = "cv")

其中~前后分别是因变量和自变量。GWR里因变量只能有1个,但自变量可以选多个,如果需要多个自变量的话,就在代码POP_SQMI之后用+号连接就行。gweight是你的空间加权的函数(随空间距离增大而不断衰减的函数,衰减率由下面要提到的带宽控制),这里用的是比较常用的高斯函数,其余的还有gwr.bisquare等函数可以调用。verbose决定是否汇报制定窗口的过程。method是决定构建带宽窗口模型的方法,这里用的cv指的是cross validation,即交叉验证法,也是最常用的方法,简单说就是把数据分成不同的组,分别用不同的方法来做回归计算,计算完了之后记录下结果,然后打乱重新分组,再回归计算,再看结果,周而复始,最后看哪种计算方法的结果最靠谱,这种方法就是最优解。还有一种很常见的选择最佳拟合模型的方法是AIC optimisation法,把method后面的cv改成aic就可以用。具体AIC optimisation是什么:AIC(赤池信息准则)_百度百科。总之,空间加权函数和带宽窗口构建方法的选择是GWR里面十分重要的步骤。

以上便是固定带宽窗口的示意图。比如我在对佐治亚做GWR,这一轮的regression target是红色的这个县,根据做出来的窗口,圆圈以内的县都要被算为红色县的邻县,其权重根据高斯函数等空间权重函数来赋值,而圆圈以外的县,空间权重都赋为0。

不喜欢固定带宽窗口的同学也可以不用它,而是用符合Tobler地理学第一定律的非固定带宽邻域统计,操作方法是在gwr.sel里面加一个命令adapt = TRUE,这样的情况下,根据你设置的k邻居数,每一轮统计的时候,和本轮对象在k以内相邻的多边形的权重参数会被赋值为0到1之间的一个数,比如下图:

我在对佐治亚做GWR,这一轮的regression target是红色的这个县,那么图上标为1的县就是红色县的1阶邻县,标为2的是2阶(邻县的邻县),标为3的是3阶(邻县的邻县的邻县)。如果用非固定带宽邻域统计,k为3,那么1、2、3都被定义为红色县的邻县,它们的权重从3到1依次增加,会按比例被赋上0和1之间的值,而其它没有标注的县,权重为0。

下一步就是用前一步做出的窗口去扫过全局区域:

gwr_result <- gwr(AVE_SIZE12 ~ POP_SQMI, data = usa_agri, bandwidth = bw,

gweight = gwr.Gauss, hatmatrix = TRUE)

这一步如果数据量大,可能会要跑一阵,跑完之后我们看看结果里面有什么:

Call:

gwr(formula = AVE_SIZE12 ~ POP_SQMI, data = usa_agri, bandwidth = bw,

gweight = gwr.Gauss, hatmatrix = TRUE)

Kernel function: gwr.Gauss

Fixed bandwidth: 205880.3

Summary of GWR coefficient estimates at data points:

Min. 1st Qu. Median 3rd Qu.Max. Global

X.Intercept. 7.3883e+01 2.1081e+02 3.2802e+02 6.6691e+02 8.5705e+03 625.5656

POP_SQMI -8.0085e+01 -4.5983e-01 -1.4704e-01 -7.3703e-02 -2.1859e-03 -0.0426

Number of data points: 3108

Effective number of parameters (residual: 2traceS - traceS'S): 119.6193

Effective degrees of freedom (residual: 2traceS - traceS'S): 2988.381

Sigma (residual: 2traceS - traceS'S): 1048.78

Effective number of parameters (model: traceS): 84.90185

Effective degrees of freedom (model: traceS): 3023.098

Sigma (model: traceS): 1042.741

Sigma (ML): 1028.4

AICc (GWR p. 61, eq 2.33p. 96, eq. 4.21): 52109.55

AIC (GWR p. 96, eq. 4.22): 52017.7

Residual sum of squares: 3287040139

Quasi-global R2: 0.4829366

基本上你做GWR该需要的结果这里都有了。比如窗口大小(Fixed bandwidth)是205880.3,意思是前一步构建的带宽窗口是半径205.88千米的圆。Effective number of parameters显示的是你带宽窗口的大小合不合适。Sigma是残差的标准差,这个值要尽量小。Residual sum of squares(RSS)也是对拟合程度的一个评估值。最重要的是最后那个R2,越靠近1说明统计的拟合度越好。我这里面Sigma很大,R2也不是很大,因为我这里只是呈现方法,用的数据本来就是互不相干、没什么太大意义的,所以不用太纠结。如果你是真正的统计数据要来做GWR,就需要注意这些值了。

然后,我们就可以把每个县的R2画在地图上。首先,前面报告里的这些数据,比如R2,要先自己去生成的GWR结果里面去找,然后自己再算一下每个县的local R2,并把它们赋值到shp里面去: