# CSVデータをウェブサイトから取得する。
d0 <- read.csv('https://stats.dip.jp/01_ds/data/load_amedas_y2022_tokyo.csv')

# データの期間を絞る（2022年5月1日～2022年6月30日）
is.term <- '2022-05-01' <= d0$date & d0$date <= '2022-06-30'
d1 <- d0[is.term, ]

n <- nrow(d1)                 # データサイズ
t <- 1:n                      # 経過時間（hour）
px <- as.POSIXct(d1$datetime) # 時間オブジェクトに変換 
y <- d1$mw / 1000             # 測定値 (単位変換：MW -> GW)
d1 <- data.frame(px = format(px, '%Y-%m-%d %H'), t, y)

head(d1) # データの先頭の5行を表示
##              px t     y
## 1 2022-05-01 00 1 21.62
## 2 2022-05-01 01 2 20.65
## 3 2022-05-01 02 3 20.61
## 4 2022-05-01 03 4 20.81
## 5 2022-05-01 04 5 20.86
## 6 2022-05-01 05 6 20.72
library('forecast')
## Registered S3 method overwritten by 'quantmod':
##   method            from
##   as.zoo.data.frame zoo
# LOESS（局所回帰）のフィッティング（fitting）
# 平滑化窓幅（span）：30日（トレンドを捉えるのに十分な窓幅を設定する）
fit1 <- loess(y ~ t,
              data = data.frame(px, t = 1:n, y), 
              span = 24 * 30 / n)

# フィット値（fitted values）がトレンドとなる。
trend <- fit1$fitted

# 折れ線プロット（原系列）
matplot(px, y, type = 'l', pch = 1, col = 2,
        main = '電力需要',
        xlab = '月 日',
        ylab = 'GW')

abline(h = 0, col = 'gray') # x軸線

# 折れ線プロット（トレンド成分）
matlines(px, 
         trend, 
         pch = 2, 
         col = 4, 
         lwd = 2)

# 凡例（はんれい）
legend('topleft', lty = 1, col = c(2, 4), 
       legend = c('原系列', 'トレンド'))

# 残余成分（トレンド除去成分）計算
detrended <- y - trend

# 折れ線プロット（残余成分）
matplot(px, detrended, type = 'l', pch = 1, col = 2,
        main = '残余成分（原系列－トレンド）',
        xlab = '月 日',
        ylab = 'GW')

abline(h = 0, col = 'gray') # x軸線

fit2 <- loess(y ~ t, 
              data = data.frame(px, t = 1:n, y = detrended), 
              span = 24 * 4 / n)

seasonal <- fit2$fitted

# 折れ線プロット（トレンド除去成分・周期成分）
matplot(px, detrended, type = 'l', pch = 1, col = 2, ylim = c(-15, 15),
        main = '残余成分（トレンド除去後の成分）・周期成分',
        xlab = '月 日',
        ylab = 'GW')

# 折れ線プロット（周期成分）
matlines(px, seasonal, pch = 2, col = 4, lwd = 2)

abline(h = 0, col = 'gray') # x軸線

# 残余成分を計算（トレンド除去成分－周期成分）
remainder <- detrended - seasonal

# 折れ線プロット（残余成分）
matplot(px, remainder, type = 'l', pch = 1, col = 2, ylim = c(-15, 15),
        main = '残余成分（トレンド，周期成分除去後の成分）',
        xlab = '月 日',
        ylab = 'GW')

abline(h = 0, col = 'gray') # x軸線

s <- mstl(x = msts(y, seasonal.period = c(24, 24*7)),
          s.window = 'periodic')

autoplot(s, main = '周期成分：時間変動なし')

s <- mstl(x = msts(y, seasonal.period = c(24, 24*7)))

autoplot(s, main = '周期成分：時間変動あり')

px <- seq(as.POSIXct('1949-01-01'),
          as.POSIXct('1960-12-01'), by = 'month')

d0 <- data.frame(px = paste(px), y = AirPassengers)

library(DT)
datatable(d0)
s <- mstl(x = msts(y, seasonal.period = c(12, 12*7)),
          s.window = 'periodic')

autoplot(s, main = '周期成分：時間変動なし')

R Markdown

This is an R Markdown document. Markdown is a simple formatting syntax for authoring HTML, PDF, and MS Word documents. For more details on using R Markdown see http://rmarkdown.rstudio.com.

When you click the Knit button a document will be generated that includes both content as well as the output of any embedded R code chunks within the document. You can embed an R code chunk like this:

summary(cars)
##      speed           dist       
##  Min.   : 4.0   Min.   :  2.00  
##  1st Qu.:12.0   1st Qu.: 26.00  
##  Median :15.0   Median : 36.00  
##  Mean   :15.4   Mean   : 42.98  
##  3rd Qu.:19.0   3rd Qu.: 56.00  
##  Max.   :25.0   Max.   :120.00

Including Plots

You can also embed plots, for example:

Note that the echo = FALSE parameter was added to the code chunk to prevent printing of the R code that generated the plot.