考古学のためのRでGIS-filter(grepl() )で必要な情報を取り出す


RでGISを行う場合に、もっとも挫折しやすいのがデータハンドリングです。手元にあるデータをそのままRで解析できるケースは非常にレアです。欲しい形に整形する作業が必ず発生します。

埋蔵文化財包蔵地の一覧表は、通常、一つの遺跡に対して複数の時期・時代が設定されています（下図）。一つのカラムに複数の時期が入力されているケースではRでの取り扱いが難しくなります。

今回は、一つのカラムに1つの時期・時代だけが入力されたデータベースの作成方法を紹介します。
![01.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/60438/c79b5240-21df-ee04-78b3-3e8b323567a8.png)


# データ読み込み

北海道教育委員会作成[「北の遺跡案内」](https://www2.wagmap.jp/hokkai_bunka/Portal)を利用して作成したシェープファイルを読み込みます。データは

```{r}
# 道南の包蔵地データ
site3100n <- st_read("isekidatabase9_JGD2000utm54.shp" ,
  options = "ENCODING=CP932") %>%  # 文字コードがShift_JISなのでエンコードを指定
  dplyr::select(sitename , period )  # 遺跡名と時代のカラムだけを取り出す
```

そのままでは、preiodカラムに複数の時期が含まれているので、時期ごとに切り分けていきます。

```{r}
>head(site3100n)
        sitename                 period                 geometry
1       虻羅遺跡 縄文(前期)、縄文(後期) POINT (405113.2 4706010)
2       最内遺跡 縄文(後期)、縄文(晩期) POINT (405571.3 4700497)
3       須築遺跡             縄文(晩期) POINT (403989.5 4716456)
4 瀬田内チャシ跡                 アイヌ   POINT (405530 4697476)
5       切梶遺跡                   不明 POINT (404885.7 4710595)
6     大光寺遺跡             縄文(中期) POINT (405770.8 4700540)
```

# `filter()`と`grepl()`を組み合わせて、時期別にデータを格納

今回の記事の肝はここです。`filter(grepl())`の組み合わせはごちゃごちゃのデータベースから情報を引き出す頼もしい関数です。
`filter(grepl( "早期" , period ) )`は、periodカラムから「早期」の文字列を含む行を取り出します。長いコードですが、同じことの繰り返しです。

```{r}
site3100n01 <-site3100n %>%
  filter(grepl( "早期" , period ) ) %>%　　# grepl()でperiodカラムに”早期”の文字列が入っている行を取り出す
  mutate( class = "早期" )  # mutate()でclassカラムを新設して、"早期"を入力
site3100n02 <-site3100n %>%
  filter(grepl( "前期" , period ) ) %>%
  mutate( class = "前期" )
site3100n03 <-site3100n %>%
  filter(grepl( "中期" , period ) ) %>%
  mutate( class = "中期" )
site3100n04 <-site3100n %>%
  filter(grepl( "後期" , period ) ) %>%
  mutate( class = "後期" )
site3100n05 <-site3100n %>%
  filter(grepl( "晩期" , period ) ) %>%
  mutate( class = "晩期" )
site3100n06 <-site3100n %>%
  filter(grepl( "続縄" , period ) ) %>%
  mutate( class = "続縄文" )
site3100n07 <-site3100n %>%
  filter(grepl( "オホーツク" , period ) ) %>%
  mutate( class = "擦文文化期" )
site3100n08 <-site3100n %>%
  filter(grepl( "擦文" , period ) ) %>%
  mutate( class = "擦文文化期" )
site3100n09 <-site3100n %>%
  filter(grepl( "中世" , period ) ) %>%
  mutate( class = "アイヌ文化期" )
site3100n10 <-site3100n %>%
  filter(grepl( "アイヌ" , period ) ) %>%
  mutate( class = "アイヌ文化期" )
site3100n11 <-site3100n %>%
  filter(grepl( "近世" , period ) ) %>%
  mutate( class = "アイヌ文化期" )
```
早期の遺跡を取り出した`site3100n01`は以下のようにclassカラムが新設され、「早期」が入力されています。

```{r}
site3100n01
         sitename                                         period
1  貝取澗山中遺跡                                     縄文(早期)
2          都遺跡 縄文(早期)、縄文(前期)、縄文(中期)、縄文(後期)
3      豊岡４遺跡                                     縄文(早期)
4      豊岡６遺跡                         縄文(早期)、縄文(中期)
5      豊岡７遺跡                                     縄文(早期)
6        松江遺跡 縄文(早期)、縄文(前期)、縄文(中期)、縄文(後期)
7        青苗遺跡   縄文(早期)、縄文(中期)、続縄文、擦文、アイヌ
8        赤川遺跡 縄文(早期)、縄文(前期)、縄文(中期)、縄文(後期)
9      砥石線遺跡                                     縄文(早期)
10       米岡遺跡                         縄文(早期)、縄文(中期)
                   geometry class
1  POINT (408201.2 4671226)  早期
2  POINT (402317.3 4675916)  早期
3  POINT (406315.6 4697112)  早期
4    POINT (406715 4696157)  早期
5  POINT (407749.6 4697158)  早期
6  POINT (373854.7 4658979)  早期
7  POINT (371605.4 4658115)  早期
8  POINT (372595.2 4658807)  早期
9  POINT (371029.8 4658900)  早期
10 POINT (371126.3 4658660)  早期
```
早期遺跡（site3100n01）をプロット

```{r}
site3100n01 %>%
  ggplot() +
    geom_sf()
```
![早期遺跡分布.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/60438/50891a5e-2dca-0fd0-4ab4-0a6b2495a408.png)

# 時期別データを`rbind()`で結合

時期ごとに切り分けたデータを`rbind()`で結合します。

```{r}
# 時期別のデータを結合
site3100mrg <- 
  rbind(site3100n01 , site3100n02 , site3100n03 , site3100n04 , site3100n05 ,
    site3100n06 , site3100n07 , site3100n08 , site3100n09 , site3100n10 ,
    site3100n11 )
```

もともとのデータの1427レコードから1990レコードに増加しています。500件ほどの複合遺跡があったということです。

```{r}
> site3100mrg
Simple feature collection with 1990 features and 3 fields
geometry type:  POINT
dimension:      XY
bbox:           xmin: 369757.2 ymin: 4584095 xmax: 514246.8 ymax: 4716456
epsg (SRID):    3100
proj4string:    +proj=utm +zone=54 +ellps=GRS80 +towgs84=0,0,0,0,0,0,0 +units=m +no_defs
```

# 時代の順序を定義して時期をfactorに変換

凡例や集計結果が時期ごとに整列するようにします。また、`st_read()`で読み込んだデータはデフォルトでcharacter（文字型）として読み込まれるので、この時点でfactorに変換します。

```{r}
# 遺跡の時期を順序を定義してfactorに変換
site3100mrg$class <- factor( site3100mrg$class ,
  levels=c( "早期" , "前期","中期","後期","晩期","続縄文","擦文文化期","アイヌ文化期") ) 
# 地図を表示
site3100mrg %>%
  ggplot() +
  geom_sf( aes( colour = class ) )
```
時期ごとに色分けされた遺跡分布図が作成できました。
![03.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/60438/34c92b2c-f3e1-01ee-8019-314b25af588c.png)

# 背景地図を設定

背景図を設定して仕上げます。分布図の場合背景図はシンプルな方が好ましく感じます。

```{r}
# OSM背景図
library(ggspatial)  # OSMを背景図にするパッケージを読み込み
site3100mrg %>%
  ggplot() +
    annotation_map_tile( type = "osm" , zoomin = 1 ) +  # OSMを背景図に読み込み　zoomin=1　が適切
    geom_sf( aes( colour = class ,) , alpha = 0.3 , show.legend = "point" ) +  # 遺跡点
    labs(caption = "\U00a9 OpenStreetMap contributors") +  # 地図の出典
    theme_void() +
    scale_colour_viridis(option = "cividis" , discrete = TRUE)
```
![04道南遺跡地図_OSM.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/60438/6a9df436-eeb8-1ac8-b46f-a2e3eaf39179.png)

```{r}
# 行政界背景図
site3100mrg %>%
  ggplot() +
    geom_sf( data = dn , fill = NA , colour = "grey70" , lwd = 0.3 ) +  # 行政界
    geom_sf( aes( colour = class ,) , alpha = 0.3 , show.legend = "point" ) +  # 遺跡点
    theme_void() +
    scale_colour_viridis( option = "cividis" , discrete = TRUE )
```
![道南遺跡地図_行政界.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/60438/98b06969-edbe-21c6-6c92-c3cb43537786.png)