空間結合と空間演算
st_join()、st_intersection()、st_buffer()、st_distance() を適用します。
「空間結合と空間演算」はCoddyKit上の無料R Academyレッスンです。 これはレッスン3/4です。 下記で完全なレッスンを無料で読むことができます。その後、ブラウザ内の組み込みコードエディタと24時間対応のAIチューターでハンズオン演習できます。 これはR Academy学習パスの一部であり、ウェブとCoddyKitアプリ全体で進捗が同期されます。 R Academyコースには全4レッスンが含まれています。
st_joinを使った空間結合
st_join(x, y)は空間的な左結合を実行します。xの各フィーチャについて、空間述語（デフォルトはst_intersects）を満たすyのフィーチャを検索し、その属性を追加します。複数のフィーチャが一致した場合は、複数の行が生成されます。
library(sf)
nc <- st_transform(
read_sf(system.file('shape/nc.shp', package = 'sf')), 4326
)
# Create random points within NC bounding box
set.seed(42)
bbox <- st_bbox(nc)
pts <- data.frame(
lon = runif(20, bbox['xmin'], bbox['xmax']),
lat = runif(20, bbox['ymin'], bbox['ymax'])
)
pts_sf <- st_as_sf(pts, coords = c('lon', 'lat'), crs = 4326)
# Spatial join: assign county attributes to each point
pts_with_county <- st_join(pts_sf, nc['NAME'])
cat('Points with county:\n')
print(head(pts_with_county[, c('NAME', 'geometry')], 5))st_within: ポイント・イン・ポリゴン
st_within(x, y)は、xの各フィーチャがyのフィーチャ内に完全に含まれるかどうかを示す論理値の疎行列を返します。st_intersectsよりも厳密であり、境界に触れているだけの場合は該当しません。
library(sf)
nc <- st_transform(
read_sf(system.file('shape/nc.shp', package = 'sf')), 4326
)
# Random points
set.seed(7)
bbox <- st_bbox(nc)
pts_sf <- st_as_sf(
data.frame(
lon = runif(30, bbox['xmin'], bbox['xmax']),
lat = runif(30, bbox['ymin'], bbox['ymax'])
),
coords = c('lon', 'lat'), crs = 4326
)
# Which points are within NC at all?
within_nc <- st_within(pts_sf, st_union(nc), sparse = FALSE)[, 1]
cat('Points within NC:', sum(within_nc), '/ 30\n')
# Keep only points inside NC
pts_inside <- pts_sf[within_nc, ]
cat('Kept:', nrow(pts_inside), 'points\n')st_intersects: 重なりの判定
st_intersects(x, y)は、フィーチャ同士が境界または内部を問わず、ジオメトリの一部を共有している場合にTRUEを返します。最も許容範囲の広いトポロジー述語です。sparse = FALSEを渡すと、インデックスベクトルのリストではなく論理行列を取得できます。
library(sf)
# Create two overlapping polygons
poly1 <- st_polygon(list(matrix(c(0,0, 2,0, 2,2, 0,2, 0,0), ncol=2, byrow=TRUE)))
poly2 <- st_polygon(list(matrix(c(1,1, 3,1, 3,3, 1,3, 1,1), ncol=2, byrow=TRUE)))
poly3 <- st_polygon(list(matrix(c(5,5, 6,5, 6,6, 5,6, 5,5), ncol=2, byrow=TRUE)))
sfc <- st_sfc(poly1, poly2, poly3, crs = 32119)
# Intersects matrix
mat <- st_intersects(sfc, sparse = FALSE)
cat('Intersects matrix:\n')
print(mat)
cat('\npoly1 intersects poly2:', mat[1, 2], '\n')
cat('poly1 intersects poly3:', mat[1, 3], '\n')st_buffer: バッファの作成
st_buffer(sf_obj, dist)は、指定した距離で各フィーチャの周囲にバッファポリゴンを作成します。地理座標系では距離の単位は度数、投影座標系ではメートルです。正確な結果を得るには、必ずバッファ作成の前に投影してください。
library(sf)
# City points in WGS84
cities <- data.frame(
city = c('Raleigh', 'Charlotte', 'Greensboro'),
lon = c(-78.639, -80.843, -79.792),
lat = c(35.779, 35.227, 36.073)
)
cities_sf <- st_as_sf(cities, coords = c('lon', 'lat'), crs = 4326)
# Project to metres for accurate buffering (NC State Plane)
cities_proj <- st_transform(cities_sf, 32119)
# 25 km buffer around each city
buffers <- st_buffer(cities_proj, dist = 25000)
cat('Buffer geometry type:', unique(as.character(st_geometry_type(buffers))), '\n')
cat('Buffer area (km2):', round(as.numeric(st_area(buffers[1, ])) / 1e6, 2), '\n')
cat('Expected (pi*25^2):', round(pi * 25^2, 2), 'km2\n')st_area: ポリゴンの面積
st_area(sf_obj)は、各ポリゴンの面積を計算します。地理座標系ではステラジアンで値を返すため、実用的ではありません。メートル単位の投影CRSでは平方メートルを返します。結果にはunitsパッケージによる単位情報が付加されます。
library(sf)
library(units)
nc <- read_sf(system.file('shape/nc.shp', package = 'sf'))
nc_proj <- st_transform(nc, 32119) # NC State Plane, metres
# Area in m2 with units
nc_proj$area_m2 <- st_area(nc_proj)
# Convert to km2
nc_proj$area_km2 <- set_units(nc_proj$area_m2, 'km^2')
# Summary statistics
cat('Smallest county:', nc_proj$NAME[which.min(nc_proj$area_km2)],
round(min(as.numeric(nc_proj$area_km2)), 1), 'km2\n')
cat('Largest county: ', nc_proj$NAME[which.max(nc_proj$area_km2)],
round(max(as.numeric(nc_proj$area_km2)), 1), 'km2\n')st_distance: ペアごとの距離
st_distance(x, y)は、xのすべてのフィーチャとyのすべてのフィーチャの間の距離行列を計算します。WGS84上では測地線距離を、投影CRS上ではユークリッド距離を使用します。結果にはunitsのメタデータが付加されます。
library(sf)
cities <- data.frame(
city = c('Raleigh', 'Charlotte', 'Greensboro', 'Wilmington'),
lon = c(-78.639, -80.843, -79.792, -77.909),
lat = c(35.779, 35.227, 36.073, 34.226)
)
cities_sf <- st_as_sf(cities, coords = c('lon', 'lat'), crs = 4326)
# Geodesic distance matrix (metres)
dist_m <- st_distance(cities_sf)
dist_km <- round(dist_m / 1000)
rownames(dist_km) <- cities$city
colnames(dist_km) <- cities$city
cat('Distance matrix (km):\n')
print(dist_km)
# Nearest neighbour for Raleigh
raleigh_dist <- as.numeric(dist_km['Raleigh', -1])
cat('\nNearest to Raleigh:',
cities$city[-1][which.min(raleigh_dist)],
min(raleigh_dist), 'km\n')st_union: ジオメトリの融合
st_union(sf_obj)は、すべてのフィーチャを単一のジオメトリに結合し、内部境界を溶解します。複数の郡ポリゴンから1つの外形ポリゴンを作成したり、複数のバッファを融合したりする場合に便利です。
library(sf)
nc <- st_transform(
read_sf(system.file('shape/nc.shp', package = 'sf')), 32119
)
# Dissolve all counties into the NC state outline
nc_state <- st_union(nc)
cat('Counties before union:', nrow(nc), '\n')
cat('Features after union:', length(nc_state), '(single geometry)\n')
cat('Geometry type:', as.character(st_geometry_type(nc_state)), '\n')
# Area check: sum of counties == state area
county_area_sum <- sum(as.numeric(st_area(nc)))
state_area <- as.numeric(st_area(nc_state))
cat('Sum of county areas:', round(county_area_sum / 1e9, 2), 'km2\n')
cat('State area: ', round(state_area / 1e9, 2), 'km2\n')st_intersection: クリッピング
st_intersection(x, y)は、各フィーチャのペアにおける幾何学的な共通部分、つまり両方のジオメトリが共有する領域を返します。一方のレイヤーをもう一方の範囲に合わせてクリッピングする際に使用します。
library(sf)
nc <- st_transform(
read_sf(system.file('shape/nc.shp', package = 'sf')), 4326
)
# Clip NC counties to a bounding rectangle (eastern NC)
clip_box <- st_as_sfc(st_bbox(c(
xmin = -78, xmax = -75,
ymin = 33.5, ymax = 36.5
), crs = 4326))
eastern_nc <- st_intersection(nc, clip_box)
cat('Original features:', nrow(nc), '\n')
cat('Clipped features: ', nrow(eastern_nc), '\n')
cat('Geometry types: ', paste(unique(as.character(st_geometry_type(eastern_nc))), collapse = ', '), '\n')dplyrを使った空間集計
st_join()をdplyr::group_by()およびsummarise()と組み合わせると、ポリゴンごとにポイントデータを集計できます。sfオブジェクトをグループ化すると、ジオメトリ列も自動的に集計されます（デフォルトでは融合されます）。
library(sf)
library(dplyr)
nc <- st_transform(
read_sf(system.file('shape/nc.shp', package = 'sf')), 4326
)
# Generate random sales events in NC
set.seed(99)
bbox <- st_bbox(nc)
pts <- st_as_sf(
data.frame(
lon = runif(200, bbox['xmin'], bbox['xmax']),
lat = runif(200, bbox['ymin'], bbox['ymax']),
sales = runif(200, 100, 1000)
),
coords = c('lon', 'lat'), crs = 4326
)
# Join then aggregate
result <- st_join(pts, nc['NAME']) |>
filter(!is.na(NAME)) |>
st_drop_geometry() |>
group_by(NAME) |>
summarise(total_sales = sum(sales), n_events = n()) |>
arrange(desc(total_sales))
cat('Top 5 counties by sales:\n')
print(head(result, 5))st_nearest_feature
st_nearest_feature(x, y)は、xの各フィーチャに対して、yで最も近いフィーチャのインデックスを返します。st_distance()と組み合わせると、最近傍フィーチャとその距離の両方を取得できます。割り当て問題に便利です。
library(sf)
# Hospitals (facilities)
hospitals <- data.frame(
name = c('Hospital A', 'Hospital B', 'Hospital C'),
lon = c(-79.5, -80.2, -78.8),
lat = c(35.8, 35.3, 36.1)
)
hosp_sf <- st_as_sf(hospitals, coords = c('lon', 'lat'), crs = 4326)
# Patient locations
patients <- data.frame(
id = 1:5,
lon = c(-79.9, -80.5, -79.0, -78.5, -80.1),
lat = c(35.7, 35.2, 35.9, 35.5, 36.0)
)
pt_sf <- st_as_sf(patients, coords = c('lon', 'lat'), crs = 4326)
# Find nearest hospital for each patient
nearest_idx <- st_nearest_feature(pt_sf, hosp_sf)
pt_sf$nearest_hospital <- hospitals$name[nearest_idx]
print(st_drop_geometry(pt_sf))実践的なワークフロー: バッファとカウント
現実の分析でよく使われるパターンは、フィーチャの集合（学校など）にバッファを作成し、その各バッファ内にいくつの関心地点（バス停など）があるかをst_join()とdplyr::count()で数える方法です。
library(sf)
library(dplyr)
# Schools as points
schools <- data.frame(
id = 1:3,
lon = c(-79.8, -80.2, -79.4),
lat = c(35.9, 35.5, 36.1)
)
schools_sf <- st_as_sf(schools, coords = c('lon', 'lat'), crs = 4326) |>
st_transform(32119) # project to metres
# Random bus stops
set.seed(5)
stops <- data.frame(
stop_id = 1:50,
lon = rnorm(50, -79.8, 0.3),
lat = rnorm(50, 35.8, 0.2)
) |>
st_as_sf(coords = c('lon', 'lat'), crs = 4326) |>
st_transform(32119)
# 2 km buffer around each school
buffers <- st_buffer(schools_sf, dist = 2000)
buffers$school_id <- schools$id
# Count bus stops within each school buffer
joined <- st_join(stops, buffers['school_id'])
counts <- count(st_drop_geometry(joined), school_id)
cat('Bus stops within 2km:\n')
print(counts)クイックチェック
レイヤーAの各ポリゴンと重なる、レイヤーBのすべてのポリゴンフィーチャを検索し、Bの属性をAに追加したいとします。どの関数を使えばよいでしょうか。
まとめ: 空間結合と空間操作
重要なポイント:
st_join(x, y)— 述語（デフォルトはst_intersects）を使った空間的な左結合st_within()— xがyの内部に完全に含まれる場合にTRUE（境界は除外）st_intersects()— xの一部がyに触れている場合にTRUE（最も許容範囲が広い）st_buffer(sf_obj, dist)— バッファを作成します。メートル単位で正確なバッファを作るには投影CRSを使用しますst_area()/st_distance()— 面積と距離を計算します（正確さのため、先に投影します）st_union()— すべてのフィーチャを融合します。st_intersection()— 一方のレイヤーをもう一方のレイヤーに合わせてクリッピングしますst_nearest_feature(x, y)— 最近傍フィーチャを検索します
library(sf)
nc <- st_transform(
read_sf(system.file('shape/nc.shp', package = 'sf')), 32119
)
# Quick operations summary
cat('Largest county:', nc$NAME[which.max(st_area(nc))], '\n')
nc_outline <- st_union(nc)
cat('State area km2:', round(as.numeric(st_area(nc_outline)) / 1e6, 0), '\n')よくある質問
「空間結合と空間演算」レッスンは無料ですか？
はい。「空間結合と空間演算」の完全なテキストはこのウェブで無料で読めます。インタラクティブに演習し（組み込みコードエディタと24時間対応のAIチューター）、R Academyコースの残りをアンロックするには、CoddyKit PROにアップグレードしてください。 R Academyコースには全4レッスンが含まれています。
「空間結合と空間演算」で何を学びますか？
st_join()、st_intersection()、st_buffer()、st_distance() を適用します。 ブラウザで直接実行するハンズオンコードでR Academyを演習し、24時間対応のAIチューターがレッスンを進める中での質問に答えます。
R Academyを始めるのに経験は必要ですか？
事前経験は必要ありません。CoddyKitのR Academyは初級者から上級者向けに構成されているため、ここから始めるか最初から始めて、自分のペースで進むことができます。 これはレッスン3/4です。
「空間結合と空間演算」レッスンにはどのくらい時間がかかりますか？
ほとんどのCoddyKitレッスンは約5～10分かかります。各レッスンはコンパクトでインタラクティブなので、着実に進歩し、ウェブとアプリ全体で正確に前回の場所から再開できます。
このR Academyレッスンでコードを書いて実行できますか？
はい。すべてのR Academyレッスンに組み込みコードエディタが含まれているため、ブラウザでリアルコードを書いて実行し、即座のAIフィードバックを取得できます。ローカル設定は不要です。