FinalReport_rmd.utf8.md

###蟲蟲危機：臺北市行道樹荔枝椿象爆發潛勢分析

###空間分析　期末報告

####森林三 B05605025 殷楷智

####森林三　B05605034　李弘恩

library(rgdal)
library(spatstat)
library(GISTools)
library(spdep)

####讀入台北市行道樹資料

setwd("D:\\Kevin Li_NTU\\107_2\\Spatial Analysis\\Final Project")
TPtree <- read.csv("./TaipeiTree_Wei.csv", encoding = "utf8")

####讀入台北市鄉鎮區

TPvill <- readOGR("./Taipei_Vill", layer = "Taipei_Vill", encoding = "utf8")

####將無患子科行道樹的資料取出並且做成SpatialPointsDataFrame

TPtree.sapind <- subset(TPtree, TPtree$Species_name_code=="DILO7"|TPtree$Species_name_code=="KOEL"|TPtree$Species_name_code=="SASA"|TPtree$Species_name_code=="LICH4")
TPtree.sapind.sp <- SpatialPointsDataFrame(coords = cbind(TPtree.sapind$x, TPtree.sapind$y), data = TPtree.sapind, proj4string = TPvill@proj4string)

####無患子科行道樹分布點圖

plot(TPvill, main = "台北市無患子科行道樹分布點圖")
points(TPtree.sapind.sp, col = "red", pch = 1, cex = 0.5)

####無患子科行道樹NNA

bb <- bbox(TPtree.sapind.sp)
x.coor <- TPtree.sapind.sp$x
y.coor <- TPtree.sapind.sp$y
x.range <- bb[1,]
y.range <- bb[2,]
sap.ppp<-ppp(x.coor,y.coor,x.range,y.range)

CIL <- envelope(sap.ppp,Lest,nsim = 99,nrank = 1)
plot(CIL,.-r~r,main = 'L funtion')

####無患子科行道樹KDE

TPvill.range <- cbind(TPvill@bbox[1,],TPvill@bbox[2,])
sap.kde <- kde2d(TPtree.sapind.sp@data$x,TPtree.sapind.sp@data$y,1000,500,TPvill.range)

image(sap.kde,asp=1,bty='n',xaxt='n',yaxt='n', main = "台北市無患子科行道樹核密度圖")
masker <- poly.outer(TPvill.range, TPvill)
add.masking(masker,col='aliceblue')
plot(TPvill, add=T)

####無患子科行道樹在各行政區的面密度

#計算行政區面積
TPvill$Area <- poly.areas(TPvill)

#算出每個行政區有多少無患子科的樹
TPvill$sapindaceae <- poly.counts(TPtree.sapind.sp, TPvill)
TPvill$sapind.den <- TPvill$sapindaceae/TPvill$Area*10^4

TPvill$VILLAGE[which.max(TPvill$sapind.den)]

#畫面量圖
lm.palette <- colorRampPalette(c("white", "orange", "red"), space = "rgb")
spplot(TPvill, zcol="sapind.den", col.regions=lm.palette(20), main="無患子科行道樹株數面積密度(棵/公頃)")

####無患子科行道樹高面密度在哪些行政區間有群聚(Cluster)，以各行政區的Gi*為檢定指標(單尾z檢定)

#找出Neighborhood
TPvill_nb <- poly2nb(TPvill, queen = T)

#轉為ListW (weighting matrix)
TPvill_nb_w <- nb2listw(TPvill_nb, zero.policy = T, style = "W")

#Standardized Gi* values
TPvill_nb_in <- include.self(TPvill_nb)
TPvill_nb_in_w <- nb2listw(TPvill_nb_in)
Gi_star <- localG(TPvill$sapind.den, TPvill_nb_in_w)

#畫圖
Gi_star_v <- as.vector(Gi_star)
quad.g <- c()
quad.g[Gi_star_v >= qnorm(0.95)] <- 1 #cluster
quad.g[Gi_star_v < qnorm(0.95)] <- 2 #non-cluster

colors.g <- c("red", "lightgrey")
plot(TPvill, border="grey", col = colors.g[quad.g], main = "Cluster Map")
legend("bottomright", legend=c("Cluster","Non-cluster"), fill=colors.g)

####有高面密度群聚現象的行政區

cbind(as.character(TPvill$TOWN[quad.g==1]), as.character(TPvill$VILLAGE[quad.g==1]))

##       [,1]     [,2]    
##  [1,] "大安區" "群英里"
##  [2,] "大安區" "虎嘯里"
##  [3,] "大安區" "臥龍里"
##  [4,] "大安區" "龍淵里"
##  [5,] "大安區" "芳和里"
##  [6,] "大安區" "敦安里"
##  [7,] "士林區" "蘭雅里"
##  [8,] "士林區" "蘭興里"
##  [9,] "士林區" "天福里"
## [10,] "士林區" "天祿里"
## [11,] "士林區" "天壽里"
## [12,] "士林區" "天和里"
## [13,] "士林區" "天山里"
## [14,] "信義區" "黎順里"
## [15,] "信義區" "黎平里"
## [16,] "大安區" "德安里"
## [17,] "大安區" "龍圖里"
## [18,] "大安區" "龍陣里"
## [19,] "大安區" "龍雲里"
## [20,] "大安區" "住安里"
## [21,] "大安區" "義安里"
## [22,] "大安區" "通安里"
## [23,] "大安區" "法治里"
## [24,] "大安區" "全安里"
## [25,] "大安區" "群賢里"
## [26,] "文山區" "興得里"
## [27,] "士林區" "義信里"
## [28,] "士林區" "前港里"
## [29,] "士林區" "岩山里"
## [30,] "士林區" "名山里"
## [31,] "士林區" "德行里"
## [32,] "士林區" "聖山里"
## [33,] "士林區" "忠誠里"
## [34,] "士林區" "芝山里"
## [35,] "士林區" "三玉里"