Mid1.utf8

109-2 空間分析　期中考一

助教杜承軒 2021.03.22

從健康促進的觀點，小學附近若有較多的速食店，容易取得高熱量食物可能影響學生的飲食健康。請以1公里範圍作為鄰近的定義，協助台北市衛生局進行以下分析。

＊初步環境建置與讀取檔案

library(sf);library(units);library(tmap);library(ggplot2);library(dplyr)
setwd("D:/1092SA/Mid1/Data/")
TP=st_read("Taipei_District.shp")
SC=st_read("Taipei_School.shp")
FF=st_read("Taipei_Fastfood.shp")
windowsFonts(TOP = windowsFont("Topedia Sans TW Beta"))

彙整列表在私立小學及公立小學1公里範圍內，加總速食店的店家數與座位數。[10%]

PV=SC[SC$type=="private",]
PB=SC[SC$type=="public",]

#way 1
PV.FF.d=st_distance(PV,FF)
PB.FF.d=st_distance(PB,FF)
FF.PV.1K=apply(PV.FF.d,2,min)<1000
FF.PB.1K=apply(PB.FF.d,2,min)<1000
A1=data.frame(type=c("private","public"),stores=c(sum(FF.PV.1K),sum(FF.PB.1K)),seats=c(sum(FF.PV.1K*FF$seat),sum(FF.PB.1K*FF$seat)))
A1

##      type stores seats
## 1 private     47  4860
## 2  public     82  8869

#way 2
PV.1K=PV%>%st_buffer(1000)%>%st_union
PB.1K=PB%>%st_buffer(1000)%>%st_union
FF.PV=FF[PV.1K,]
FF.PB=FF[PB.1K,]
A1=data.frame(type=c("private","public"),stores=c(nrow(FF.PV),nrow(FF.PB)),seats=c(sum(FF.PV$seat),sum(FF.PB$seat)))
A1

##      type stores seats
## 1 private     47  4860
## 2  public     82  8869

將各行政區內速食店的總座位數，除以該行政區小學學生總人數，可得出「每位學生平均可分配到速食店座位數」，該數值越大表示學生對於速食店的可及性越大。請繪製主題地圖，以面量圖(choropleth map)表示該數值在台北市各行政區的分布（需繪製面量圖圖例）。[10%]

#way 1
TP.SC=st_contains(TP,SC,sparse=F) # = st_is_within_distance(TP,SC,0)
TP.FF=st_contains(TP,FF,sparse=F)
TP$student=TP.SC%*%SC$student %>% as.numeric
TP$seat=TP.FF%*%FF$seat %>% as.numeric
TP$access=TP$seat/TP$student
A2=tm_shape(TP) + tm_polygons("access",title="可近性指標")+tm_layout(frame = F, fontfamily = "TOP")
A2

#way 2
SC.TP=st_intersection(SC,TP)
FF.TP=st_intersection(FF,TP)
TPxStudent=xtabs(student~Code,SC.TP)
TPxStudent=data.frame(Code=as.numeric(names(TPxStudent)),Student=as.numeric(TPxStudent))
TP=left_join(TP,TPxStudent,by = "Code")
TPxSeat=xtabs(seat~Code,FF.TP)
TPxSeat=data.frame(Code=as.numeric(names(TPxSeat)),Seat=as.numeric(TPxSeat))
TP=left_join(TP,TPxSeat,by = "Code")
TP$Access=TP$Seat/TP$Student
A2=qtm(TP, fill="Access",fill.title="可近性指標")+tm_layout(frame = F, fontfamily = "TOP")
A2

繪製泡泡地圖(Bubble Map)，標示1公里範圍內有最多間速食店之前八名的小學分布位置，並調整圓圈的大小，來表現這些學校1公里內的速食店家座位數總和（需繪製圓圈大小的圖例）。(20%)[10%]

#way 1
SC.FF=st_is_within_distance(SC,FF,1000)
SC$store=lapply(SC.FF,length)%>%unlist
SC$seat=lapply(SC.FF,function(x) sum(FF$seat[x]))%>%unlist

#way 2
SC.FF.1K=st_is_within_distance(SC,FF,1000,sparse=F)
SC$store=rowSums(SC.FF.1K)
SC$seat=SC.FF.1K%*%FF$seat # = apply(SC.FF.1K,1,function(x) sum(x*FF$seat))

SC8=SC[order(SC$store,decreasing=T)[1:8],]
SC8$鄰近速食店座位總數=SC8$seat
A3=tm_shape(TP)+tm_borders()+tm_shape(SC8)+tm_bubbles("鄰近速食店座位總數",col="blue")+tm_layout(frame = F, fontfamily = "TOP")
A3

#way 3
SC.buf1K=st_buffer(SC,1000)
SC1K.FF=st_intersection(SC.buf1K,FF)
SC.FF.x=xtabs(cbind(1,seat)~SID,SC1K.FF)
SC.FF.df=as.data.frame.matrix(SC.FF.x)
colnames(SC.FF.df)=c("Store","Seat")
SC.FF.df$SID=SC.FF.df%>%row.names%>%as.numeric
SC=left_join(SC,SC.FF.df)

SC8=SC[order(SC$store,decreasing=T)[1:8],]
SC8$附近速食店座位總數=SC8$seat
A3=tm_shape(TP)+tm_borders()+tm_shape(SC8)+tm_bubbles("附近速食店座位總數",col="blue")+tm_layout(frame = F, fontfamily = "TOP")
A3

#補充：which.max
x=SC$store
ids=c()
for(i in 1:8){
  id=which.max(x)
  x[id]=NA
  ids=c(ids,id)
}
all(ids==order(SC$store,decreasing=T)[1:8])

## [1] TRUE

衛生局想列出優先輔導飲食健康的行政區。
(1) 繪製盒狀圖(box plot)，評估各行政區內在學校附近的速食店之座位數。[20%]
X軸：行政區
Y軸：該行政區內距離學校1公里內的速食店之座位數

#way 1
FF$District=TP$District[apply(TP.FF,2,which)]
FF.1K=FF[apply(SC.FF.1K,2,any),]
A4.1=ggplot(FF.1K)+geom_boxplot(aes(x=District,y=seat,fill=District))+xlab("行政區")+ylab("速食店座位數")+theme_minimal()+theme(text=element_text(family="TOP"),legend.position = "none")
A4.1

#way 2
FFTP.SC1K=st_within(FF.TP,SC.buf1K) # = st_is_within_distance(FF,SC.buf1K,0)
FFTP.1K=FF.TP[lapply(FFTP.SC1K,length)!=0,]
A4.1=ggplot(FF.1K)+geom_boxplot(aes(x=District,y=seat,fill=District))+xlab("行政區")+ylab("速食店座位數")+theme_minimal()+theme(text=element_text(family="TOP"),legend.position = "none")
A4.1

(2) 用合適的統計檢定方法，比較各行政區平均值之間是否有統計上的顯著差異。需列出假設檢定與結論。[10%]

設定顯著水準為\({\alpha=0.05}\)
\(H_0:\)各行政區鄰近學校之速食店座位數的平均值相同
\(H_a:\)各行政區鄰近學校之速食店座位數的平均值有異
檢定多母體平均值是否相同，使用AVOVA檢定。

A4.2=aov(seat~District,FF.1K)%>%summary
A4.2

##             Df Sum Sq Mean Sq F value   Pr(>F)    
## District    11  86115    7829    7.87 3.76e-09 ***
## Residuals   82  81567     995                     
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

\(F=7.87\)，\(p-value=0\)。
\(p-value<\alpha\)，拒絕虛無假設。
結論：各行政區鄰近學校之速食店座位數的平均值有異。

設定0.5公里內；0.5～1公里；1～1.5公里；1.5～2公里四種不同的環域範圍，並將台北市的小學分成北區與南區（如下）。
北區：北投、士林、大同、中山、松山、內湖
南區：萬華、中正、大安、信義、南港、文山
(1) 繪製折線圖(line chart)，比較北區與南區的小學，不同環域範圍內的速食店店家總數。[20%]

SC$District=TP$District[apply(TP.SC,2,which)]
#或
SC$dist=c()
for(i in SC$SID) SC$dist[i]=SC.TP$District[SC.TP$SID==i]

SC$Region="S"
SC$Region[SC$District%in%c("北投區","士林區","大同區","中山區","松山區","內湖區")]="N"
#SC$Region[SC$District%in%c("萬華區","中正區","大安區","信義區","南港區","文山區")]="S"
N.SC=SC[SC$Region=="N",]
S.SC=SC[SC$Region=="S",]

#way1
N.FF=st_distance(N.SC,FF)
N.Count=table(cut(apply(N.FF,2,min),c(0:4*500,Inf)))[1:4]
S.FF=st_distance(S.SC,FF)
S.Count=table(cut(apply(S.FF,2,min),c(0:4*500,Inf)))[1:4]
NS.Count=data.frame(region=rep(c("N","S"),each=4),dist=rep(seq(.5,2,.5),2),Count=c(N.Count,S.Count))
A5.1=ggplot(NS.Count)+geom_line(aes(x=dist,y=Count,colour=region),size=1)+geom_point(aes(x=dist,y=Count,colour=region),size=3)+scale_color_manual("台北市分區",values=c(2,4),labels = c("北區","南區"))+xlim(0,2)+ylim(0,50)+xlab("環域距離")+ylab("搜尋區內速食店家數")+theme_minimal()+theme(text=element_text(family="TOP"),legend.position = "bottom")
A5.1

#way2
FF.buf=function(pts,dist=seq(.5,2,.5)){
  count=c()
  for(d in dist){
  buf=pts%>%st_buffer(set_units(d,km))%>%st_union
  count=c(count,nrow(FF[buf,]))
  }
  count=count-c(0,count[-length(count)])
  return(count)
}

N.count=FF.buf(N.SC);S.count=FF.buf(S.SC)
NS.count=data.frame(region=rep(c("N","S"),each=4),dist=rep(seq(.5,2,.5),2),count=c(N.count,S.count))
A5.1=ggplot(NS.count)+geom_line(aes(x=dist,y=count,colour=region),size=1)+geom_point(aes(x=dist,y=count,colour=region),size=3)+scale_color_manual("台北市分區",values=c(2,4),labels = c("北區","南區"))+xlim(0,2)+ylim(0,50)+xlab("環域距離")+ylab("搜尋區內速食店家數")+theme_minimal()+theme(text=element_text(family="TOP"),legend.position = "bottom")
A5.1

(2) 用合適的統計檢定方法，比較北區與南區的小學，在不同環域範圍內的速食店店家總數是否有一致的趨勢。需列出假設檢定與結論。[10%]

設定顯著水準為\({\alpha=0.05}\)
\(H_0:\)北區和南區在不同環域範圍的趨勢相同
\(H_a:\)北區和南區在不同環域範圍的趨勢不同
檢定兩組計數資料的分布狀況，使用卡方檢定。

A5.2=chisq.test(cbind(N.count,S.count))
A5.2

## 
##  Pearson's Chi-squared test
## 
## data:  cbind(N.count, S.count)
## X-squared = 2.6062, df = 3, p-value = 0.4564

\(\chi^2=2.606\)，\(p-value=0.456\)。
\(p-value>\alpha\)，不拒絕虛無假設。
結論：北區和南區在不同環域範圍的趨勢相同。

109-2 空間分析 期中考一

助教 杜承軒 2021.03.22

109-2 空間分析　期中考一

助教杜承軒 2021.03.22