GGRS: Geoscience, GIS, and Remote Sensing

지구과학, GIS, 그리고 원격탐사 블로그입니다.

REMOTE SENSING

PyTorch: TorchGeo를 이용한 지리공간 딥러닝 (2)

유병혁 2023. 1. 22. 23:29

안녕하세요? 이번 글은 시리즈 글의 두번째로 'PyTorch의 TorchGeo를 이용한 지리공간 딥러닝' 실습을 다뤄보겠습니다. 이 글은 지리공간/금융 데이터 사이언티스트 Maurício Cordeiro(마우리시오 코르데이로) 님이 작성하신 'PyTorch의 TorchGeo를 이용한 지리공간 분석용 인공지능(2부)' 내용을 정리 및 보완한 것입니다. 원본 글은 다음 링크를 참고하시면 됩니다. 좋은 글을 공유해주신 마우리시오 코르데이로 님께 감사의 마음을 전합니다.

[1] Artificial Intelligence for Geospatial Analysis with Pytorch’s TorchGeo (Part 1)
[2] Artificial Intelligence for Geospatial Analysis with Pytorch’s TorchGeo (Part 2)
[3] Artificial Intelligence for Geospatial Analysis with Pytorch’s TorchGeo (Part 3)

 

Artificial Intelligence for Geospatial Analysis with Pytorch’s TorchGeo (Part 1)

An end-to-end deep learning geospatial segmentation project using Pytorch and TorchGeo packages

towardsdatascience.com

제가 작성한 첫번째 글은 다음 링크를 참고하시면 됩니다.

 

PyTorch: TorchGeo를 이용한 지리공간 딥러닝 (1)

안녕하세요? 이번 글은 시리즈 글의 첫번째로 'PyTorch의 TorchGeo를 이용한 지리공간 딥러닝' 실습을 다뤄보겠습니다. 이 글은 지리공간/금융 데이터 사이언티스트 Maurício Cordeiro(마우리시오 코르데

foss4g.tistory.com

첫번째 글에서 TorchGeo를 사용하여 RasterDataset을 생성하고 RandomSampler를 사용하여 이미지 샘플을 추출하는 방법을 학습했었습니다. 이번에는 패치와 해당 레이블을 튜플로 다루기 위해 이미지와 마스크를 IntersectionDataset으로 연결해 보겠습니다.

 

환경 설정

Colab에서 환경 설정을 위해 다음 패키지(packaging, rasterio, torchgeo)들을 추가 설치합니다.

# 패키지 추가 설치
!pip install packaging -q
!pip install rasterio -q
!pip install torchgeo -q

Colab의 내장 패키지 matplotlib을 최신 버전으로 업그레이드합니다. 관련해서 다음 링크 글을 참고하시면 좋겠습니다.

 

Colab: 내장 패키지 확인 및 업그레이드 방법 소개

안녕하세요? 이번에는 Google Colab에서 내장 패키지 목록을 확인하고, 패키지 버전을 업그레이드 하는 방법을 간략히 소개해 보겠습니다. 이번 글은 인도네시아 PT Mobility Doctor Indonesia에서 개발자

foss4g.tistory.com

# 내장 패키지 업그레이드
from packaging import version

import matplotlib
print("matplotlib 현재 버전: " + matplotlib.__version__)

if(version.parse(matplotlib.__version__) < version.parse("3.6.0")):
    print("matplotlib 업그레이드")
    !pip install --upgrade matplotlib -q
    import os
    os.kill(os.getpid(), 9) # 런타임 종료

Google 드라이브에 위치한 Earth Surface Water Dataset(dset-s2.zip)을 액세스 후 압축 해제합니다. 코드 중 assert문은 claim과 같은 뜻을 가지고 있습니다. 'dset-s2' 경로가 없으면 False 값이 반환되어 AssertionError가 발생합니다.

# Google 드라이브 엑세스 후 데이터셋 압축 해제
from google.colab import drive
from pathlib import Path

drive.mount('/gdrive', force_remount=True)
!unzip -qq "/gdrive/My Drive/GEODATA/dset-s2.zip"
root = Path('dset-s2')
assert root.exists()

RasterDataset 클래스

 

훈련 이미지, 훈련 마스크, 검증 이미지, 검증 마스크에 대한 RasterDataset 객체를 생성합니다. 사용하는 데이터셋의 좌표계는 EPSG:3395로 일치하게 하고 해상도는 10m, 이미지 값은 DN(화소)를 10,000으로 나눈 Reflectance(반사율)로 변환합니다. 여기서 중요한 것은 훈련 마스크, 검증 마스크는 'image'가 아니라는 것을 RasterDataset에 알리는 것입니다. 이렇게 하면 샘플이 마스크 데이터셋에서 추출될 때 일반적으로 이미지에 사용되는 "image" 대신 "mask" 키가 결합된 데이터셋으로 반환됩니다. 훈련 이미지는 64장, 검증 이미지는 31장입니다.

# 패키지 불러오기
import rasterio as rio
import logging

logger = logging.getLogger("rasterio")
logger.setLevel(logging.ERROR) # ERROR 미만 로깅 메시지 무시

import torchgeo
from torchgeo.datasets import RasterDataset, unbind_samples, stack_samples
from torchgeo.samplers import RandomGeoSampler, Units
from torch.utils.data import DataLoader

# Reflectance(반사율) = DN(화소) / 10000
def scale(item: dict):
    item['image'] = item['image'] / 10000
    return item

train_imgs = RasterDataset(root=(root/'tra_scene').as_posix(), crs='epsg:3395', res=10, transforms=scale) # 훈련 이미지(반사율)
train_msks = RasterDataset(root=(root/'tra_truth').as_posix(), crs='epsg:3395', res=10) # 훈련 마스크(레이블)

valid_imgs = RasterDataset(root=(root/'val_scene').as_posix(), crs='epsg:3395', res=10, transforms=scale) # 검증 이미지(반사율)
valid_msks = RasterDataset(root=(root/'val_truth').as_posix(), crs='epsg:3395', res=10) # 검증 마스크(레이블)

train_msks.is_image = False # 훈련 마스크 이미지 여부: x
valid_msks.is_image = False # 검증 마스크 이미지 여부: x

print(train_imgs, valid_imgs) # 훈련 이미지, 검증 이미지(64장, 31장)
output:
RasterDataset Dataset
    type: GeoDataset
    bbox: BoundingBox(minx=-12324108.040448364, maxx=18180258.35959181, miny=-3799980.7640827284, maxy=11631943.161442367, mint=0.0, maxt=9.223372036854776e+18)
    size: 64 RasterDataset Dataset
    type: GeoDataset
    bbox: BoundingBox(minx=-8570422.010074684, maxx=16205540.179445159, miny=-4640708.937046933, maxy=9086073.349200383, mint=0.0, maxt=9.223372036854776e+18)
    size: 31

 

Sampler 클래스

샘플러는 데이터셋으로부터 훈련에 필요한 더 작은 고정 크기의 패치를 만드는 역할입니다. 여기서는 TorchGeo에서 제공하는 샘플러 중에서 RandomGeoSampler 클래스를 사용합니다. 샘플러는 원본 이미지에서 고정 크기의 임의 경계상자(bbox)를 선택하고, 이 bbox를 사용하여 데이터셋 중 원하는 이미지 부분을 질의합니다.

 

예시로 살펴볼까요?! 훈련 이미지 데이터셋에서 512*512 크기 샘플 1장을 추출합니다. 경계상자를 사용하여 훈련 이미지 중 해당 부분을 질의할 수 있습니다. 샘플 유형은 [6, 512, 512]입니다. 여기서 6은 밴드 개수를 뜻합니다.

sampler = RandomGeoSampler(train_imgs, size=512, length=1) # 512*512 크기로 1장 추출
bbox = next(iter(sampler)) # 경계상자
sample = train_imgs[bbox] # 샘플
sample['image'].shape # 유형
output:
torch.Size([6, 512, 512])

샘플을 플로팅해본 결과입니다.

import torch
import matplotlib.pyplot as plt
plt.rcParams["figure.figsize"] = (7,7)

arr = torch.clamp(sample['image'], min=0, max=1).numpy()
rgb = arr.transpose(2, 1, 0)[:, :, [2, 1, 0]] # Natural Color

plt.imshow(rgb*3)

이런 방식으로 훈련 이미지 샘플은 130장, 검증 이미지 샘플은 64장을 추출합니다. 데이터셋은 더 작은 크기인 512*512 패치로 잘려질 것입니다.

train_sampler = RandomGeoSampler(train_imgs, size=512, length=130, units=Units.PIXELS) # 훈련 이미지 샘플
valid_sampler = RandomGeoSampler(valid_imgs, size=512, length=64, units=Units.PIXELS)  # 검증 이미지 샘플

 

IntersectionDataset 클래스

이미지와 마스크를 튜플로 다루기 위해 IntersectionDataset을 생성합니다. 다음과 같이 & 연산자를 사용하여 간단하게 이미지와 마스크(레이블)를 결합할 수 있습니다.

#IntersectionDataset 클래스
train_dset = train_imgs & train_msks # 훈련 데이터셋(이미지/마스크) 결합
valid_dset = valid_imgs & valid_msks # 검증 데이터셋(이미지/마스크) 결합

 

DataLoader 클래스

DataLoader는 샘플 집단, 즉 배치를 신경망 훈련 절차에 제공하는 역할을 합니다. 여기서 batch_size는 각 배치의 샘플 수(=8장)를, collate_fn은 여러 샘플을 하나의 단일 배치로 연결하는 방법(=stack_samples)을 지정합니다. DataLoader는 반복하여 배치를 가져올 수 있습니다. 실험용으로 첫번째 배치를 불러와 봅니다. 

# 데이터 로더
train_dataloader = DataLoader(train_dset, sampler=train_sampler, batch_size=8, collate_fn=stack_samples)
valid_dataloader = DataLoader(valid_dset, sampler=valid_sampler, batch_size=8, collate_fn=stack_samples)

train_batch = next(iter(train_dataloader)) # 단일 훈련 배치 (샘플 8장)
valid_batch = next(iter(valid_dataloader)) # 단일 검증 배치 (샘플 8장)
print(train_batch.keys(), valid_batch.keys())
output:
dict_keys(['image', 'crs', 'bbox', 'mask']) dict_keys(['image', 'crs', 'bbox', 'mask'])

 

 

Batch Visualization

훈련 배치를 가시화해 보겠습니다. 단일 훈련 배치는 훈련 이미지 8장, 훈련 마스크 8장으로 구성되어 있습니다. 배치 플롯 결과는 다음과 같습니다.

from typing import Iterable, List

# 이미지 플롯
def plot_imgs(images: Iterable, axs: Iterable, chnls: List[int] = [2, 1, 0], bright: float = 3.):
    for img, ax in zip(images, axs):
        arr = torch.clamp(bright * img, min=0, max=1).numpy()
        rgb = arr.transpose(1, 2, 0)[:, :, chnls]
        ax.imshow(rgb)
        ax.axis('off')

# 마스크 플롯
def plot_msks(masks: Iterable, axs: Iterable):
    for mask, ax in zip(masks, axs):
        ax.imshow(mask.squeeze().numpy(), cmap='Blues')
        ax.axis('off')

# 이미지/마스크 배치 플롯
def plot_batch(batch: dict, bright: float = 3., cols: int = 4, width: int = 5, chnls: List[int] = [2, 1, 0]):
    samples = unbind_samples(batch.copy()) # 단일 배치를 개별 샘플로 전환
    n = 2 * len(samples) if ('image' in batch) and ('mask' in batch) else len(samples)
    rows = n//cols + (1 if n%cols != 0 else 0)
    _, axs = plt.subplots(rows, cols, figsize=(cols*width, rows*width))

    # 격자 생성
    if ('image' in batch) and ('mask' in batch):
        # 짝수 축: 이미지 플롯
        plot_imgs(images=map(lambda x: x['image'], samples), axs=axs.reshape(-1)[::2], chnls=chnls, bright=bright)
        # 홀수 축: 마스크 플롯
        plot_msks(masks=map(lambda x: x['mask'], samples), axs=axs.reshape(-1)[1::2])
    else:
        if 'image' in batch:
            plot_imgs(images=map(lambda x: x['image'], samples), axs=axs.reshape(-1), chnls=chnls, bright=bright)   
        elif 'mask' in batch:
            plot_msks(masks=map(lambda x: x['mask'], samples), axs=axs.reshape(-1))

# 이미지/마스크 배치 플롯(각8장)
plot_batch(train_batch)

 

Data Normlization(Standardization)

일반적으로 머신러닝(딥러닝 포함)은 피처 스케일링(feature scaling)의 이점을 얻습니다. 여기서는 데이터를 데이터의 평균으로 뺀 값을 표준편차로 나누어 표준화(X' = (X - 평균) / (표준편차)) 해보겠습니다. 이렇게 하면 평균은 0이고 표준편차가 1인 표준 정규분포를 따릅니다. 일단 계산을 위해 데이터셋 각 6개 채널에 대한 평균과 표준편차를 구합니다. 데이터셋을 메모리에 한번에 적재하는 부담을 피하고자, 단일 이미지에 대한 평균과 표준편차를 반복 누적한 후 이미지 개수로 나누어 반환하는 함수를 정의합니다.

# 전체 데이터셋에 대한 통계(평균, 표준편차) 계산
def calc_statistics(dset: RasterDataset):
        # RasterDataset을 파일 리스트로 변환
        # 예: files[0] = 'dset-s2/tra_scene/S2A_L2A_20190125_N0211_R034_6Bands_S2.tif'
        files = [item.object for item in dset.index.intersection(dset.index.bounds, objects=True)]

        # 초기화(평균, 표준편차)
        accum_mean = 0
        accum_std = 0

        # 단일 이미지 통계 누적
        for file in files:
            img = rio.open(file).read()/10000
            accum_mean += img.reshape((img.shape[0], -1)).mean(axis=1)
            accum_std += img.reshape((img.shape[0], -1)).std(axis=1)

        # 이미지 개수로 나누어 반환
        return accum_mean / len(files), accum_std / len(files)

훈련 이미지 6개 채널에 대한 평균과 표준편차는 다음과 같습니다.

mean, std = calc_statistics(train_imgs)
print(mean, std)
output:
[0.0771449  0.09890421 0.09758993 0.22216185 0.1854808  0.13288888] [0.04496952 0.05038998 0.06053346 0.10840577 0.0993342  0.08219175]

평균과 표준편차를 구했으니 해당 값을 이용해 데이터를 정규화하는 클래스를 정의합니다.

forward 메소드는 정규화를 수행합니다. 아래 그림과 같이 플로팅을 위해 정규화를 되돌리는 revert 메소드도 정의합니다. 

normalize = MyNormalize(*calc_statistics(train_imgs))
norm_batch = normalize.forward(train_batch) # 정규화
plot_batch(norm_batch)

revert는 다음과 같이 확인할 수 있습니다.

batch = normalize.revert(norm_batch) # 정규화 취소
plot_batch(batch)

 

nn.Sequential 클래스

TorchGeo는 신경망의 성능 개선을 위해 위성 이미지에서 분광지수를 쉽게 생성하고 배치에 추가할 수 있습니다. 각 배치에 분광지수와 정규화를 추가합니다. 여기서는 배치에 센티널 2호 위성 이미지로부터 계산한 NDWI(Normalized Difference Water Index, 정규수분지수), MNDWI(Modified NDWI, 보정수분지수), NDVI(Normalized Difference Vegetation Index, 정규식생지수)를 추가합니다.

 

변환된 배치는 9개 채널(3종 분광지수, 정규화된 6개 채널)이 있습니다. 분광지수는 정규화된 값을 사용하지 않기 때문에 normalize가 가장 마지막 줄에 배치되었습니다.

from torchgeo.transforms import indices

tfms = torch.nn.Sequential(
    indices.AppendNDWI(index_green=1, index_nir=3),  # NDWI(Normalized Difference Water Index, 정규수분지수)
    indices.AppendNDWI(index_green=1 , index_nir=5), # MNDWI(Modified Normalized Difference. Water Index, 보정수분지수)
    indices.AppendNDVI(index_nir=3, index_red=2),    # NDVI(Normalized Difference Vegetation Index, 정규식생지수)
    normalize # 정규화
)

transformed_batch = tfms(train_batch)
print(transformed_batch['image'].shape, transformed_batch['mask'].shape)
output:
torch.Size([8, 9, 512, 512]) torch.Size([8, 1, 512, 512])

이번 글은 원본 데이터셋의 이미지와 마스크를 결합하여 IntersectionDataset을 생성하는 법, nn.Sequential 클래스를 사용하여 원본 데이터에 변환을 추가하는 법을 학습해 봤습니다. 다음 글에서는 훈련 루프, 손실 함수를 만들고 새로 만든 딥러닝 결과를 확인해 보겠습니다. 

 

PyTorch-TorchGeo를_이용한_지리공간 딥러닝_2.ipynb

Colaboratory notebook

colab.research.google.com