R语言空间分析：学校选址，实现代码（链接有完整工程），Rstudio

R语言空间分析完整工程

期末作业，要求可以参考如下：（环境：RStudio），请用RStudio打开，直接用R可能乱码。 R语言实现代码如下：

library(readr) library(sp) library(maptools) library(rgdal) library(rgeos) library(ggplot2) library(GISTools) library(raster) library(recharts) library(REmap) library(leaflet) library(pdftables) library(ggthemes) library(reshape2) setwd("C:/Users/96233/Desktop/R语言期末作业/期末工程") #任务一：爱尔兰小学教育资源配置的时空发变化 #-------------------------------------------------------------------------------------------------------------- #读入数据 这一板块的代码无需运行，生成的文件已经保存在文件夹里了， #爱尔兰学校资源统计年鉴数据 这个API Key只能有50页免费翻译，我快用完了，好像还有十来页，你可以跑一个转换看看效果。 library(pdftables) #从PDF读入，将pdf转为csv。(将5个pdf都转换为csv，有免费次数限制) #2010-2011 write.csv(head(iris, 20), file = "Output/2010_2011.csv", row.names = FALSE) get_remaining("djfvsm9il2li")#这个是在https://pdftables.com注册申请的API convert_pdf("Data/schools-Key-Statistics-2010-2011.pdf","Output/2010_2011.csv",api_key = "djfvsm9il2li") #2012 write.csv(head(iris, 20), file = "Output/2012.csv", row.names = FALSE) get_remaining("djfvsm9il2li")#这个是在https://pdftables.com注册申请的API convert_pdf("Data/schools-Key-Statistics-2012.pdf","Output/2012.csv",api_key = "djfvsm9il2li") #2012-2013 write.csv(head(iris, 20), file = "Output/2012-2013.csv", row.names = FALSE) get_remaining("djfvsm9il2li")#这个是在https://pdftables.com注册申请的API convert_pdf("Data/schools-Key-Statistics-2012-2013.pdf","Output/2012-2013.csv",api_key = "djfvsm9il2li") #2013-2014 write.csv(head(iris, 20), file = "Output/2013-2014.csv", row.names = FALSE) get_remaining("djfvsm9il2li")#这个是在https://pdftables.com注册申请的API convert_pdf("Data/schools-Key-Statistics-2013-2014.pdf","Output/2013-2014.csv",api_key = "djfvsm9il2li") #2014-2015 write.csv(head(iris, 20), file = "Output/2014-2015.csv", row.names = FALSE) get_remaining("djfvsm9il2li")#这个是在https://pdftables.com注册申请的API convert_pdf("Data/schools-Key-Statistics-2014-2015.pdf","Output/2014-2015.csv",api_key = "djfvsm9il2li") #-------------------------------------------------------------------------------------------------------------- #由于年鉴数据较多，我们选择其中的某些部分做分析 #1、总体教育资源变化 School_stat_2010_2011=read.csv("Output/总体情况/2010-2011.csv") School_stat_2011_2012=read.csv("Output/总体情况/2012.csv") School_stat_2012_2013=read.csv("Output/总体情况/2012-2013.csv") School_stat_2013_2014=read.csv("Output/总体情况/2013-2014.csv") School_stat_2014_2015=read.csv("Output/总体情况/2014-2015.csv") #三个等级教育人数变化 first_level<-c(School_stat_2010_2011[1,2],School_stat_2011_2012[1,2],School_stat_2012_2013[1,3],School_stat_2013_2014[1,3],School_stat_2014_2015[1,3]) second_level<-c(School_stat_2010_2011[4,2],School_stat_2011_2012[4,2],School_stat_2012_2013[4,3],School_stat_2013_2014[4,3],School_stat_2014_2015[4,3]) third_level<-c(School_stat_2010_2011[8,2],School_stat_2011_2012[8,2],School_stat_2012_2013[12,3],School_stat_2013_2014[12,3],School_stat_2014_2015[12,3]) first_level=as.numeric(first_level) second_level=as.numeric(second_level) third_level=as.numeric(third_level) #定义坐标轴 time<-c("2010","2011","2012","2013","2014") #画图 plot(1:5,first_level,pch=50,type='l',col="red",xaxt="n",xlab="Year",ylab="Number",ylim=c(100000,900000)) axis(1,labels = time,at=1:5) lines(1:5,second_level,pch=50,type='l',col="green") lines(1:5,third_level,pch=50,type='l',col="blue") legend("topright",legend=c("first_level","second_level","third_level"),col=c("red","green","blue"),lty=1,lwd=2) #2、区域资源变化 Area_stat_2010_2011=read.csv("Output/地区分布/2010-2011.csv") Area_stat_2011_2012=read.csv("Output/地区分布/2012.csv") Area_stat_2012_2013=read.csv("Output/地区分布/2012-2013.csv") Area_stat_2013_2014=read.csv("Output/地区分布/2013-2014.csv") Area_stat_2014_2015=read.csv("Output/地区分布/2014-2015.csv") #2010-2011 County<-c() Pupils<-c() Classes<-c() Area_stat_2010_2011$County<-as.character(Area_stat_2010_2011$County) #使用饼状图表现空间 for(i in 1:33) { County<-c(County,Area_stat_2010_2011["County"][i,1]) Pupils<-c(Pupils,Area_stat_2010_2011["Pupils"][i,1]) Classes<-c(Classes,Area_stat_2010_2011["Classes"][i,1]) } rm(i) pie(Pupils,labels = County,main="2010-2011爱尔兰学生分布情况（各郡）",radius = 1,cex.mian=2,cex=0.6) pie(Classes,labels = County,main="2010-2011爱尔兰班级分布情况（各郡）",radius = 1,cex.mian=2,cex=0.6) #2011-2012 County<-c() Pupils<-c() Classes<-c() Area_stat_2011_2012$County<-as.character(Area_stat_2011_2012$County) #使用饼状图表现空间 for(i in 1:33) { County<-c(County,Area_stat_2011_2012["County"][i,1]) Pupils<-c(Pupils,Area_stat_2011_2012["Pupils"][i,1]) Classes<-c(Classes,Area_stat_2011_2012["Classes"][i,1]) } rm(i) pie(Pupils,labels = County,main="2011-2012爱尔兰学生分布情况（各郡）",radius = 1,cex.mian=2,cex=0.6)#学生人数 pie(Classes,labels = County,main="2011-2012爱尔兰班级分布情况（各郡）",radius = 1,cex.mian=2,cex=0.6)#班级数目 #后面几年原理一样，复制粘贴改一下名字就可以 #-------------------------------------------------------------------------------------------------------------- #任务二：制作2002-2006年人口变化图，找出人口增长最大的几个ED #1、研究区域内ED尺度下的人口信息 DublinEDs<-readShapeSpatial('Data/DublinEDs.shp') #访问其中data：DublinEDs@data$(列名) ED_name=as.character(DublinEDs@data$DED_NAME) #行政区名字 ED_population2002=as.numeric(DublinEDs@data$POP02)#2002年人口 ED_population2006=as.numeric(DublinEDs@data$POP06)#2006年人口 pop_data<-data.frame(ED_name,ED_population2002,ED_population2006) library(ggthemes) library(reshape2) #包含melt() melt_pop_data<-melt(pop_data,id.vars = "ED_name",variable.name = "Year",value.name = "Number") ggplot(melt_pop_data,aes(ED_name,Number,fill=Year))+geom_bar(stat="identity",position="dodge")+theme(axis.title.x = element_text(size=1,angle = 0)) #人口变化图得到，记得把宽度扩大200倍以能够看到清楚的横坐标 #2、找出人口增长最大的几个ED ED_population_change2002_2006<-c() for(i in 1:322) { ED_population_change2002_2006<-c(ED_population_change2002_2006,ED_population2006[i]-ED_population2002[i]) } #读出最大的10个ED的序号 index_out=order(ED_population_change2002_2006,decreasing = TRUE)[1:10] rm(i) EDs_add_most=ED_name[index_out] EDs_add_most #-------------------------------------------------------------------------------------------------------------- #任务三：交通可达性地图 #1、研究区域的主要道路网络 DublinRoads<-readShapeSpatial("Data/DublinRoads.shp",verbose = TRUE) plot(DublinRoads,pch=20,col="red") #生成缓冲区 DublinRoads.buf=gBuffer(DublinRoads,width=3000) plot(DublinRoads.buf,col="green",add=T) #2、研究区域的ED DublinEDs<-readShapeSpatial('Data/DublinEDs.shp') plot(DublinEDs,pch=20,add=T) #3、Dublin市区小学位置与描述 Dublin_PrimarySchools_all=read.csv("Data/Dublin_PrimarySchools.csv") #生成空间点 Spt_schools<-SpatialPoints(list(Dublin_PrimarySchools_all['EASTING_'],Dublin_PrimarySchools_all['NORTHING'])) plot(Spt_schools,col="blue",add=T) plot(DublinRoads,pch=20,col="red",add=T) legend("topleft",legend=c("Boundary","Road","School","3km Buffer"),col=c("black","red","blue","green"),lty=1,lwd=2,cex = 1) title("交通可达性地图") #-------------------------------------------------------------------------------------------------------------- #任务四：选择3个新学校的位置分布图 注意plot的顺序，会按图层顺序覆盖的：先画缓冲区，再画EDs，再画路 #要求：1、新学校距离主路不超过3km：建立主路3km缓冲区。 #生成缓冲区 DublinRoads<-readShapeSpatial("Data/DublinRoads.shp",verbose = TRUE) DublinRoads.buf=gBuffer(DublinRoads,width=3000) plot(DublinRoads.buf,col="green") #区域 DublinEDs<-readShapeSpatial('Data/DublinEDs.shp',verbose = TRUE) plot(DublinEDs,pch=20,add=T) #画路 plot(DublinRoads,pch=20,col="red",add=T) #由于画图的关系，先把缓冲区所确定的区域求出来。 ED_intersect=intersect(DublinRoads.buf,DublinEDs) plot(ED_intersect,col="maroon") plot(DublinRoads,pch=20,col="yellow",add=T)#重新画路 #要求：2、新学校区域内没有多宗教的小学 #读入该区域小学数据 Dublin_PrimarySchools_all=read.csv("Data/Dublin_PrimarySchools.csv") Dublin_PrimarySchools_multi=Dublin_PrimarySchools_all[Dublin_PrimarySchools_all$DENOMINATI=="MULTI DENOMINATIONAL",] #生成多宗教学校的空间点位置 Spt_schools_multi<-SpatialPoints(list(Dublin_PrimarySchools_multi['EASTING_'],Dublin_PrimarySchools_multi['NORTHING'])) plot(Spt_schools_multi,col="blue",add=T) #读取多宗教小学所在的区域 flag=FALSE#判断点是否在区域中，0表示不在区域内 Selected_EDs_index<-c() for(i in 1:322) { for(j in 1:22) { flag=gContains(DublinEDs[i,2],Spt_schools_multi[j]) #此处调试的时候发现，这样可以表示出符合gcontains参数的sp类型的polygon if(flag==TRUE) { Selected_EDs_index<-c(Selected_EDs_index,i) } } } rm(flag) rm(i) rm(j) #Selected_EDs_index里面存的就是有multi的区域的序号 for(i in 1:22) { polygon(DublinEDs@polygons[[Selected_EDs_index[i]]]@Polygons[[1]]@coords,col="white") } rm(i) plot(DublinRoads,pch=20,col="green",add=T)#重新画路 #要求：3、新学校位置所在的区域内的已有小学的招生人数应当相对较少，使资源利用最大化。 #由于不知道每个学校在哪个区，所以用循环判断每个学校在那个区域ED #所有学校的空间点 Spt_schools_all<-SpatialPoints(list(Dublin_PrimarySchools_all['EASTING_'],Dublin_PrimarySchools_all['NORTHING'])) EDs_school_index<-c()#保存每一所学校所在的ED的索引 flag=FALSE for(i in 1:421) { for(j in 1:322) { flag=gContains(DublinEDs[j,2],Spt_schools_all[i]) if(flag==TRUE) { EDs_school_index<-c(EDs_school_index,j) j=1 break } } } rm(i) rm(j) #取得学校对应的ED索引后，为了方便归类，做一个data.frame EDs_school_pupils_num=Dublin_PrimarySchools_all["ENROLLMENT"] ED_school_dataframe=data.frame(EDs_school_index,EDs_school_pupils_num)#EDs_school_index是原始表中的小学所在的ED在ED表中的索引序号 #归类结果保存至info_list info_list=aggregate(x=ED_school_dataframe$ENROLLMENT,by=list(ED_school_dataframe$EDs_school_index),FUN=sum) info_school_num=info_list[1:200,2]#为了方便排序而提取 info_school_num_index_out=order(info_school_num,decreasing = TRUE)[1:190]#找出招生较多的190个，此向量里的索引对应info_list #此处的思想是，由于图层操作的原因，不太好直接将招生人数较少的ED标出，所以我选择倒着 #选出最多的190个，在图层上做减法，会方便些。 for(i in 1:190) { j=EDs_school_index[info_school_num_index_out[i]]#由于EDs_school_index中存的是小学对应的ED的索引序号，故直接在ED表中寻找。 #让我们看看效果 polygon(DublinEDs@polygons[[j]]@Polygons[[1]]@coords,col="white") } rm(i) rm(j) #此时剩下的着色区域代表：1、缓冲区与ED的交集；2、除去有多宗教小学的区域； #3、除去教育资源最多的190个区域（总共有200个ED有小学数据），等效于凸显出教育资源最少的10个区域和无教育资源的区域 #剩下的工作，可以人为在筛选出来的区域中均匀选取3个，进行新学校的位置选择。

完整项目工程压缩包地址。