からだにいいもの Rで解析：メモリと時間を節約できる「cheapr」パッケージの紹介
Rで大規模なデータを扱う際、処理速度の低下やメモリの不足が分析の妨げになることは少なくありません。特に標準的な基本関数をそのまま使うと、データ量が増えるにつれて計算コストが課題となる場面があります。cheaprパッケージは、こうした課題を解決するために、標準コマンドよりも高速でメモリ効率の良いコマンドを多数提供します。手軽にコードを置き換えるだけで、より効率的なプログラミングを実現できるのが大きな強みです。
パッケージバージョンは1.5.1。Windows 11 x64 (build 26200)のR version 4.6.1で確認しています。
＜おすすめのRに関する書籍です＞
Rグラフィックスクックブック 第2版 ―ggplot2によるグラフ作成のレシピ集 | Winston Chang, 石井弓美子, 河内崇
ggplot2でのグラフ作成をレシピ形式でまとめた一冊。この図をどう描くか、という引き方ができる。
パッケージのインストール
下記コマンドを実行してください。
# パッケージのインストール
install.packages("cheapr")
# パッケージの読み込み
library("cheapr")
# コマンド例実行で必要なパッケージを読み込み
library("bench")スポンサーリンク
コマンド例
詳細はコメント、パッケージのヘルプを確認してください。
as_discrete() の使用例
cut()の代替関数です。数値を区間ごとにグループ化する点は同じですが、区間の書式（開始点・終了点の表示形式や、区切り記号）を関数として柔軟に指定できるのが特徴です。年齢層や得点レンジなど、「見やすいラベルの階級分け」を作りたいときに重宝します。
as_discrete() のオプション（全12個）:
| 引数 | 意味 | 初期値 |
|---|---|---|
| x | 数値ベクトル | |
| breaks | 区切り点（ブレークポイント） | |
| left_closed | 左側閉区間にするか、右側閉区間にするか | |
| include_endpoint | 端点を含めるかどうか | FALSE |
| include_oob | 範囲外の値を含めるかどうか | FALSE |
| ordered | 結果を順序付き因子にするかどうか | FALSE |
| intv_start_fun | 区間の開始点の書式を整えるための関数 | |
| intv_end_fun | 区間の終了点の書式を整えるための関数 | |
| intv_closers | 閉区間であることを表す記号を指定する、長さ2の文字ベクトル | |
| intv_openers | 開区間であることを表す記号を指定する、長さ2の文字ベクトル | |
| intv_sep | 区間の開始点と終了点を区切るための、長さ1の文字ベクトル | |
| inf_label | 無限大を含む区間に使うラベル |
# as_discrete()はcut()とよく似ているが、
# 離散区間ごとに書式関数や記号を指定できる分、より柔軟
# 書式関数を使って年齢層をきれいに分類する例
ages <- 1:100
age_group <- function(x, breaks){
age_groups <- as_discrete(
x,
breaks = breaks,
intv_sep = "-",
intv_end_fun = function(x) x - 1,
intv_openers = c("", ""),
intv_closers = c("", ""),
include_oob = TRUE,
ordered = TRUE
)
# 以下は最後の年齢層のラベルを付け替えているだけ
lvls <- levels(age_groups)
n_lvls <- length(lvls)
max_ages <- paste0(max(breaks), "+")
attr(age_groups, "levels") <- c(lvls[-n_lvls], max_ages)
age_groups
}
age_group(ages, seq(0, 80, 20))
age_group(ages, seq(0, 25, 5))
age_group(ages, 5)
# as_discrete()でcut()の挙動をほぼ再現するには次のようにする
cheapr_cut <- function(x, breaks, right = TRUE,
include.lowest = FALSE,
ordered.result = FALSE){
if (length(breaks) == 1){
breaks <- get_breaks(x, breaks, pretty = FALSE,
expand_min = FALSE, expand_max = FALSE)
adj <- diff(range(breaks)) * 0.001
breaks[1] <- breaks[1] - adj
breaks[length(breaks)] <- breaks[length(breaks)] + adj
}
as_discrete(x, breaks, left_closed = !right,
include_endpoint = include.lowest,
ordered = ordered.result,
intv_start_fun = function(x) formatC(x, digits = 3, width = 1),
intv_end_fun = function(x) formatC(x, digits = 3, width = 1))
}
x <- rnorm(100)
cheapr_cut(x, 10)
identical(cut(x, 10), cheapr_cut(x, 10))実行結果：
[1] 0-19 0-19 0-19 0-19 0-19 0-19 0-19 0-19 0-19 0-19 0-19 0-19
[13] 0-19 0-19 0-19 0-19 0-19 0-19 0-19 20-39 20-39 20-39 20-39 20-39
[25] 20-39 20-39 20-39 20-39 20-39 20-39 20-39 20-39 20-39 20-39 20-39 20-39
[37] 20-39 20-39 20-39 40-59 40-59 40-59 40-59 40-59 40-59 40-59 40-59 40-59
[49] 40-59 40-59 40-59 40-59 40-59 40-59 40-59 40-59 40-59 40-59 40-59 60-79
[61] 60-79 60-79 60-79 60-79 60-79 60-79 60-79 60-79 60-79 60-79 60-79 60-79
[73] 60-79 60-79 60-79 60-79 60-79 60-79 60-79 80+ 80+ 80+ 80+ 80+
[85] 80+ 80+ 80+ 80+ 80+ 80+ 80+ 80+ 80+ 80+ 80+ 80+
[97] 80+ 80+ 80+ 80+
Levels: 0-19 < 20-39 < 40-59 < 60-79 < 80+
[1] 0-4 0-4 0-4 0-4 5-9 5-9 5-9 5-9 5-9 10-14 10-14 10-14
[13] 10-14 10-14 15-19 15-19 15-19 15-19 15-19 20-24 20-24 20-24 20-24 20-24
[25] 25+ 25+ 25+ 25+ 25+ 25+ 25+ 25+ 25+ 25+ 25+ 25+
[37] 25+ 25+ 25+ 25+ 25+ 25+ 25+ 25+ 25+ 25+ 25+ 25+
[49] 25+ 25+ 25+ 25+ 25+ 25+ 25+ 25+ 25+ 25+ 25+ 25+
[61] 25+ 25+ 25+ 25+ 25+ 25+ 25+ 25+ 25+ 25+ 25+ 25+
[73] 25+ 25+ 25+ 25+ 25+ 25+ 25+ 25+
…(以下略)gcd() の使用例
最大公約数（gcd）・最小公倍数（scm）を求める関数です。Rには専用の関数が用意されておらず、自作するか、他パッケージに頼ることが多い処理ですが、cheaprならベクトル全体に対して高速に計算できます。さらに2つの数値専用のバイナリ法版（gcd2()／scm2()）も用意されており、こちらはユークリッドの互除法よりもさらに高速に動作します。
gcd() のオプション（全6個）:
| 引数 | 意味 | 初期値 |
|---|---|---|
| x | 数値ベクトル | |
| tol | 許容誤差 | sqrt(.Machine$double.eps) |
| na_rm | TRUEの場合、NA値を無視する | TRUE |
| round | TRUEの場合、出力をround(gcd, digits)の形で丸める。digitsはtolの値から自動的に決まる | TRUE |
| break_early | 実験的な機能で、浮動小数点数にのみ適用される。TRUEの場合、計算を早期に打ち切る | TRUE |
| y | 数値ベクトル |
# バイナリ版
# 15と25の最大公約数をバイナリ法で計算
gcd2(15, 25)
# 15と5から25までの5刻みの数列の最大公約数をバイナリ法で計算
gcd2(15, seq(5, 25, 5))
# 15と5から25までの5刻みの数列の最小公倍数をバイナリ法で計算
scm2(15, seq(5, 25, 5))
# 15と25の最小公倍数をバイナリ法で計算
scm2(15, 25)
# ベクトル全体の最大公約数
# ベクトル(0, 5, 25)の最大公約数を計算
gcd(c(0, 5, 25))
# ベクトル(0, 5, 25)の最大公約数の結果をマーク
mark(gcd(c(0, 5, 25)))
# 10^5個の正規分布に従う乱数を生成してxに代入
x <- rnorm(10^5)
# 乱数ベクトルxの最大公約数を計算
gcd(x)
# 乱数ベクトルxの最大公約数を丸め処理なしで計算
gcd(x, round = FALSE)
# 乱数ベクトルxの最大公約数の結果をマーク
mark(gcd(x))実行結果：
[1] 5
[1] 5 5 15 5 5
[1] 15 30 15 60 75
[1] 75
[1] 5
# A tibble: 1 × 13
expression min median `itr/sec` mem_alloc `gc/sec` n_itr n_gc total_time
<bch:expr> <bch> <bch:> <dbl> <bch:byt> <dbl> <int> <dbl> <bch:tm>
1 gcd(c(0, 5, … 700ns 800ns 1158668. 0B 0 10000 0 8.63ms
# ℹ 4 more variables: result <list>, memory <list>, time <list>, gc <list>
[1] 1.5e-08
[1] 1.490116e-08
# A tibble: 1 × 13
expression min median `itr/sec` mem_alloc `gc/sec` n_itr n_gc total_time
<bch:expr> <bch:tm> <bch:> <dbl> <bch:byt> <dbl> <int> <dbl> <bch:tm>
1 gcd(x) 700ns 800ns 1146828. 0B 0 10000 0 8.72ms
# ℹ 4 more variables: result <list>, memory <list>, time <list>, gc <list>c_() の使用例
c()の高速版です。特に、リストにまとめた大量のオブジェクトを結合する場面（do.call(c, list)に相当する処理）で威力を発揮します。オプションの.args引数を使って直接リストを渡すと、do.call()経由で呼び出すよりも効率的に結合できます。
c_() のオプション（全1個）:
| 引数 | 意味 | 初期値 |
|---|---|---|
| .args | …（可変長引数）の代わりに使う引数 | NULL |
# c()と同じように結合する
c_(1, 2, 3:5)
# 初期設定では列ではなく行を結合する
c_(new_df(x = 1:3), new_df(x = 4:10))
# 結合したいオブジェクトのリストがある場合は
# do.callより効率的な.argsを使う
list_of_objs <- rep_(list(0), 10^4)
bench::mark(
do.call(c, list_of_objs),
do.call(c_, list_of_objs),
c_(.args = list_of_objs) # 最速
)実行結果：
[1] 1 2 3 4 5
x
1 1
2 2
3 3
4 4
5 5
6 6
7 7
8 8
9 9
10 10
# A tibble: 3 × 13
expression min median `itr/sec` mem_alloc `gc/sec` n_itr n_gc total_time
<bch:expr> <bch:t> <bch:> <dbl> <bch:byt> <dbl> <int> <dbl> <bch:tm>
1 do.call(c,… 194.3µs 229µs 4083. 78.2KB 184. 1507 68 369ms
2 do.call(c_… 778.9µs 968µs 961. 156.3KB 129. 358 48 372ms
3 c_(.args =… 96.2µs 106µs 9003. 78.2KB 10.8 4171 5 463ms
# ℹ 4 more variables: result <list>, memory <list>, time <list>, gc <list>興味深いのは、同じc_()でもdo.call()を経由するとdo.call()自体のオーバーヘッドが乗ってしまい、base Rのdo.call(c, list)よりむしろ遅くなる点です。.args引数に直接リストを渡す書き方が、速度・メモリの両面で最も効率的でした。
＜おすすめのRに関する書籍です＞
Rで楽しむベイズ統計入門［しくみから理解するベイズ推定の基礎］ Data Science Library
ベイズ推定のしくみを、仕組みの理解から丁寧に追う入門書。
www.amazon.co.jp
overview() の使用例
summary()とstr()を合わせて拡張したような関数です。数値変数・カテゴリ変数それぞれの要約統計量を一括表示でき、hist = TRUEにするとミニヒストグラム（インラインのスパークライン）も併せて表示されるため、データの分布をひと目で確認できます。
overview() のオプション（全3個）:
| 引数 | 意味 | 初期値 |
|---|---|---|
| x | ベクトルまたはデータフレーム | |
| digits | 要約統計量を表示する際の小数点以下の桁数 | getOption(“cheapr.digits”, 2) |
| hist | インラインヒストグラムを表示するかどうか | FALSE |
# 要約統計量を表示
overview(iris)
# ヒストグラム付き
overview(airquality, hist = TRUE)
# 小数点以下0桁に丸める
overview(airquality, digits = 0)
# overview()全体に適用されるオプションを設定できる
options(cheapr.digits = 1)
overview(rnorm(100))
options(cheapr.digits = 2) # 初期値に戻す実行結果：
obs: 150
cols: 5
----- Numeric -----
col n_missng p_complt n_unique mean p0 p25 p50 p75 p100 iqr sd
1 Sepal.Length 0 1 35 5.84 4.3 5.1 5.8 6.4 7.9 1.3 0.83
2 Sepal.Width 0 1 23 3.06 2 2.8 3 3.3 4.4 0.5 0.44
3 Petal.Length 0 1 43 3.76 1 1.6 4.35 5.1 6.9 3.5 1.77
4 Petal.Width 0 1 22 1.2 0.1 0.3 1.3 1.8 2.5 1.5 0.76
hist
1 ▆▇▇▅▂
2 ▁▆▇▂▁
3 ▇▁▆▇▂
4 ▇▁▇▅▃
----- Categorical -----
col n_missng p_complt n_unique n_levels min max
1 Species 0 1 3 3 setosa virginica
obs: 153
cols: 6
…(以下略)copy() の使用例
Rの「コピー・オン・モディファイ」の挙動を、明示的に制御するための関数群です。shallow_copy()は属性を複製せずに参照だけをコピーする軽量版、deep_copy()は属性も含めて完全に複製するフル版です。大きなオブジェクトを関数内で書き換えたいが元データは壊したくない、という場面で使い分けます。
copy() のオプション（全1個）:
| 引数 | 意味 | 初期値 |
|---|---|---|
| x | コピーの対象となるオブジェクト |
df <- new_df(x = sample.int(10^4))
# メモリ割り当てに注目
mark(shallow_copy(df), iterations = 1)
mark(deep_copy(df), iterations = 1)
# どちらの場合もdfのアドレスは変化する
address(df);address(shallow_copy(df));address(deep_copy(df))
# シャローコピーでは属性は複製されない
address(attr(df, "names"));address(attr(shallow_copy(df), "names"))
# ディープコピーでは属性も複製される
address(attr(df, "names"));address(attr(deep_copy(df), "names"))
# シャローコピーの有無によるインプレースでの属性追加の違い
invisible(attrs_add(df, key = TRUE, .set = TRUE))
attr(df, "key")
# インプレースで属性を削除
invisible(attrs_add(df, key = NULL, .set = TRUE))
# シャローコピーの場合
invisible(attrs_add(shallow_copy(df), key = TRUE, .set = TRUE))
# 'key'属性はシャローコピー側にのみ追加され、元のdfには追加されない
attr(df, "key")実行結果：
# A tibble: 1 × 13
expression min median `itr/sec` mem_alloc `gc/sec` n_itr n_gc total_time
<bch:expr> <bch> <bch:> <dbl> <bch:byt> <dbl> <int> <dbl> <bch:tm>
1 shallow_copy… 2.4µs 2.4µs 416667. 2.25KB 0 1 0 2.4µs
# ℹ 4 more variables: result <list>, memory <list>, time <list>, gc <list>
# A tibble: 1 × 13
expression min median `itr/sec` mem_alloc `gc/sec` n_itr n_gc total_time
<bch:expr> <bch:> <bch:> <dbl> <bch:byt> <dbl> <int> <dbl> <bch:tm>
1 deep_copy(d… 23.7µs 23.7µs 42194. 41.3KB 0 1 0 23.7µs
# ℹ 4 more variables: result <list>, memory <list>, time <list>, gc <list>
[1] "0000018d22404820"
[1] "0000018d23107428"
[1] "0000018d230ae238"
[1] "0000018d22404ac0"
[1] "0000018d22404ac0"
[1] "0000018d22404ac0"
[1] "0000018d230b29e0"
[1] TRUE
NULLfactor_() の使用例
factor()の高速・高機能版です。レベルの並べ替え（頻度順など）や、上位n件だけを残して残りを「その他」にまとめるlevels_lump()など、カテゴリカルデータの前処理でよく使う操作がひとまとまりになっています。
factor_() のオプション（全14個）:
| 引数 | 意味 | 初期値 |
|---|---|---|
| x | 対象となるベクトルを指定 | integer() |
| levels | オプションの因子レベル | NULL |
| order | 因子レベルをソートするかどうか | TRUE |
| na_exclude | NA値を因子レベルから除外するかどうか | TRUE |
| ordered | 結果を順序付き因子にするか | is.ordered(x) |
| where | NAレベルを配置する位置 | |
| name | NAレベルの名称を指定 | |
| order_by | order_byの中央値に基づく並べ替え用のベクトル | |
| decreasing | 並べ替え後のレベルを降順にするかどうか | |
| .fun | 各レベルに適用するリネーム用の関数 | |
| n | 計算対象とする上位n個のレベル数 | |
| prop | 計算対象とする上位の割合 | |
| other_category | 「その他」カテゴリーの名称 | |
| ties | 同順位（タイ）の処理方法 |
x <- factor_(sample(letters[sample.int(26, 10)], 100, TRUE), levels = letters)
x
# 使用中／未使用のレベル
levels_used(x)
levels_unused(x)
# 未使用のレベルを削除
levels_drop(x)
# 出現頻度上位3文字
lumped_letters <- levels_lump(x, 3)
levels_count(lumped_letters)
# 「その他」カテゴリーを削除するにはlevels_rm()を使う
levels_count(levels_rm(lumped_letters, "Other"))
# levels_lump()を使えば、因子でないデータに対しても
# 汎用的な上位n件抽出関数を作れる
get_top_n <- function(x, n){
f <- levels_lump(factor_(x, order = FALSE), n = n)
levels_count(f)
}
get_top_n(x, 3)
# 因子のレベルを出現頻度順に並べる方法は次のとおり：
levels(levels_lump(x, prop = 1)) # 多い順
levels(levels_lump(x, prop = -1)) # 少ない順実行結果：
[1] z c e k c t t e c e x x y x k k k t w k a a x r z a e y y e z z w c z x t
[38] e t k a x z r e w w t w a t w e r z r t a y r e e e y x t x c r r k y e w
[75] w z z z t z y y r a k z y c y w z a k z x w k k a k
Levels: a b c d e f g h i j k l m n o p q r s t u v w x y z
[1] "a" "c" "e" "k" "r" "t" "w" "x" "y" "z"
[1] "b" "d" "f" "g" "h" "i" "j" "l" "m" "n" "o" "p" "q" "s" "u" "v"
[1] z c e k c t t e c e x x y x k k k t w k a a x r z a e y y e z z w c z x t
[38] e t k a x z r e w w t w a t w e r z r t a y r e e e y x t x c r r k y e w
[75] w z z z t z y y r a k z y c y w z a k z x w k k a k
Levels: a c e k r t w x y z
name count prop
1 z 14 0.14
2 e 12 0.12
3 k 12 0.12
4 Other 62 0.62
name count prop
1 z 14 0.3684211
2 e 12 0.3157895
3 k 12 0.3157895
name count prop
…(以下略)
is_na() の使用例
is.na()やcolSums(is.na(x))などを組み合わせて自作していたNA集計処理を、num_na()・row_na_counts()・col_na_counts()といった専用関数でまとめて高速に行えます。行・列単位での欠損値カウントや、欠損だらけの行の抽出・除去が簡単に書けます。
is_na() のオプション（全3個）:
| 引数 | 意味 | 初期値 |
|---|---|---|
| x | ベクトル、リスト、データフレームまたは行列 | |
| recursive | リストに対して再帰的に処理を適用するかどうか | FALSE |
| names | 行名・列名を追加するかどうか |
### データの準備########
x <- 1:10
x[c(1, 5, 10)] <- NA
num_na(x)
which_na(x)
which_not_na(x)
row_nas <- row_na_counts(airquality, names = TRUE)
col_nas <- col_na_counts(airquality, names = TRUE)
row_nas
col_nas
df <- sset(airquality, j = 1:2)
########
# データ内のNAの数
num_na(df)
# どの行が空か
row_na <- row_all_na(df)
sset(df, row_na)
# 空の行を削除
sset(df, which_(row_na, invert = TRUE))
# または
na_rm(df)
# または
sset(df, row_na_counts(df) < ncol(df))実行結果：
[1] 3
[1] 1 5 10
[1] 2 3 4 6 7 8 9
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20
0 0 0 0 2 1 0 0 0 1 1 0 0 0 0 0 0 0 0 0
21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40
0 0 0 0 1 1 2 0 0 0 0 1 1 1 1 1 1 0 1 0
41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
0 1 1 0 1 1 0 0 0 0 0 1 1 1 1 1 1 1 1 1
61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80
1 0 0 0 1 0 0 0 0 0 0 1 0 0 1 0 0 0 0 0
81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100
0 0 1 1 0 0 0 0 0 0 0 0 0 0 0 1 1 1 0 0
101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120
0 1 1 0 0 0 1 0 0 0 0 0 0 0 1 0 0 0 1 0
121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
141 142 143 144 145 146 147 148 149 150 151
…(以下略)sequences() の使用例
seq()やsequence()の拡張版です。開始点や増分をベクトルで複数指定して一括生成できるほか、移動平均などのローリング集計で使う「ウィンドウ／ラグの添字列」を作るユーティリティ（window_sequence()・lag_sequence()・lead_sequence()）も含まれています。
sequences() のオプション（全9個）:
| 引数 | 意味 | 初期値 |
|---|---|---|
| size | 数列の長さを指定するベクトル | |
| from | 数列の開始値 | |
| by | 数列の増分（公差） | |
| add_id | 数列のID番号を名前として追加するかどうか | FALSE |
| as_list | 数列をリストとして返すかどうか。TRUEにすると各数列が個別の要素として返る | |
| to | 数列の終了値 | |
| k | ウィンドウ／ラグのサイズ | |
| partial | 部分的なウィンドウ／ラグを返すかどうか | TRUE |
| ascending | ウィンドウの数列を昇順にするかどうか | TRUE |
# この2つの関数は似た働きをする
sequence(1:3);sequence_(1:3)
# sequence_()は任意の数値ベクトルの数列に対応できる
sequence(1:3, by = 0.1);sequence_(1:3, by = 0.1)
# 数列のリストとして返すこともできる
sequence_(1:3, by = 0.1, as_list = TRUE)
# 数列にIDを追加する
sequence_(1:3, by = 0.1, add_id = TRUE)
# これを手早くデータフレームに変換する
seqs <- sequence_(1:3, by = 0.1, add_id = TRUE)
new_df(name = names(seqs), seq = seqs)
sequence(c(3, 2), by = c(-0.1, 0.1));sequence_(c(3, 2), by = c(-0.1, 0.1))
# seq()のベクトル化版
seq_(1, 10, by = c(1, 0.5))
# 上記と同じ結果
c(seq(1, 10, 1), seq(1, 10, 0.5))
# 同様に、数列のリストとして
# 2つの異なる開始点と2つの異なる増分を指定
seq_(from = c(-1, 1), 3, by = c(1, 0.5), as_list = TRUE)
# 最終的な数列の長さを求めるにはseq_sizeコマンド
sizes <- seq_size(1, 10, by = c(1, 0.5))
print(paste(c("sequence sizes: (", sizes, ") total size:", sum(sizes)),
collapse = " "))
# 数列のリストとして返す
# この長さは上記のコードの結果と一致する
seq_(1, 10, by = c(1, 0.5), as_list = TRUE) |>
list_lengths()
# 異なる増分での日付の数列
from <- Sys.Date()
to <- from + 10
by <- c(1, 2, 3)
date_seqs <- seq_(from, to, by, as_list = TRUE)
lapply(date_seqs, function(x) `class<-`(x, "Date"))
# ローリング計算
# 長さ10のベクトルに対するサイズ3のウィンドウ数列
# 過去方向を見たときのウィンドウサイズを示す
window_sequence(10, 3, partial = FALSE)
window_sequence(10, 3, partial = TRUE)
window_sequence(c(3, 5), 3)
window_sequence(c(3, 5), 3, partial = FALSE)
window_sequence(c(3, 5), 3, partial = TRUE, ascending = FALSE)
# 長さ10のベクトルに対するサイズ3のラグ数列
lag_sequence(10, 3, partial = FALSE)
# 未来方向にどれだけ見るか
lead_sequence(10, 3, partial = FALSE)
lag_sequence(10, 3, partial = TRUE)
lead_sequence(10, 3, partial = TRUE)実行結果：
[1] 1 1 2 1 2 3
[1] 1 1 2 1 2 3
[1] 1 1 1 1 1 1
[1] 1.0 1.0 1.1 1.0 1.1 1.2
[[1]]
[1] 1
[[2]]
[1] 1.0 1.1
[[3]]
[1] 1.0 1.1 1.2
1 2 2 3 3 3
1.0 1.0 1.1 1.0 1.1 1.2
name seq
1 1 1.0
2 2 1.0
3 2 1.1
4 3 1.0
…(以下略)base Rとの速度比較まとめ
記事内で紹介した関数について、base Rの相当処理と実際に速度を比較してみました。benchパッケージの代わりに、下記の自作の簡易ベンチマーク関数（proc.time()ベースで、処理をまとめて実行してから1回あたりの時間を算出する方式）を使っています。
# 簡易ベンチマーク関数（時計の分解能が粗い環境でも
# 安定した計測値が得られるよう、複数回をまとめて実行してから
# 1回あたりの時間を算出する）
simple_mark <- function(..., reps = 20) {
exprs <- as.list(substitute(list(...)))[-1]
labels <- sapply(exprs, deparse)
envir <- parent.frame()
fmt_time <- function(sec) {
if (sec < 1e-6) sprintf("%.0fns", sec * 1e9)
else if (sec < 1e-3) sprintf("%.2f\u00b5s", sec * 1e6)
else if (sec < 1) sprintf("%.2fms", sec * 1e3)
else sprintf("%.2fs", sec)
}
results <- lapply(exprs, function(e) {
invisible(eval(e, envir)) # ウォームアップ
# 1回の所要時間からバッチサイズを自動決定
t0 <- proc.time()[["elapsed"]]
invisible(eval(e, envir))
one_time <- proc.time()[["elapsed"]] - t0
batch <- if (one_time < 1e-4) 2000L
else if (one_time < 1e-3) 500L
else if (one_time < 1e-2) 50L
else if (one_time < 1e-1) 5L
else 1L
per_call <- numeric(reps)
for (i in seq_len(reps)) {
t0 <- proc.time()[["elapsed"]]
for (j in seq_len(batch)) invisible(eval(e, envir))
per_call[i] <- (proc.time()[["elapsed"]] - t0) / batch
}
per_call
})
data.frame(
expression = labels,
min = sapply(results, function(x) fmt_time(min(x))),
median = sapply(results, function(x) fmt_time(median(x)))
)
}
# ① as_discrete() vs cut()（10万件の乱数）
# ※ cheapr_cut()はas_discrete()のセクションで定義した自作関数。
cheapr_cut <- function(x, breaks, right = TRUE,
include.lowest = FALSE,
ordered.result = FALSE){
if (length(breaks) == 1){
breaks <- get_breaks(x, breaks, pretty = FALSE,
expand_min = FALSE, expand_max = FALSE)
adj <- diff(range(breaks)) * 0.001
breaks[1] <- breaks[1] - adj
breaks[length(breaks)] <- breaks[length(breaks)] + adj
}
as_discrete(x, breaks, left_closed = !right,
include_endpoint = include.lowest,
ordered = ordered.result,
intv_start_fun = function(x) formatC(x, digits = 3, width = 1),
intv_end_fun = function(x) formatC(x, digits = 3, width = 1))
}
set.seed(1)
x <- rnorm(10^5)
simple_mark(cut(x, 10), cheapr_cut(x, 10), reps = 20)
# ② factor_() vs factor()（10万件の文字ベクトル）
set.seed(1)
y <- sample(letters, 10^5, replace = TRUE)
simple_mark(factor(y), factor_(y), reps = 20)
# ③ col_na_counts() vs base Rでの列ごとのNA集計
set.seed(1)
df <- as.data.frame(matrix(rnorm(10^5), ncol = 100))
df[sample(length(df[[1]]), 500), 1] <- NA
df[sample(length(df[[2]]), 300), 2] <- NA
base_colsums <- function(d) colSums(is.na(d))
simple_mark(base_colsums(df), col_na_counts(df), reps = 20)
# ④ sequences() vs base Rでの複数開始点の数列一括生成（200回）
base_seqs <- function() lapply(1:200, function(i) seq(i, i + 20, by = 0.5))
cheapr_seqs <- function() seq_(1:200, 1:200 + 20, by = 0.5, as_list = TRUE)
simple_mark(base_seqs(), cheapr_seqs(), reps = 20)
# ⑤ c_()の呼び出し方の違い（1万要素のリスト結合）
list_of_objs <- rep_(list(0), 10^4)
simple_mark(
do.call(c, list_of_objs),
do.call(c_, list_of_objs),
c_(.args = list_of_objs),
reps = 15
)
# ⑥ shallow_copy() vs deep_copy()
df2 <- new_df(x = sample.int(10^4))
simple_mark(shallow_copy(df2), deep_copy(df2), reps = 30)実行結果：
simple_mark(cut(x, 10), cheapr_cut(x, 10), reps = 20)
expression min median
1 cut(x, 10) 2.25ms 2.34ms
2 cheapr_cut(x, 10) 1.26ms 1.30ms
simple_mark(factor(y), factor_(y), reps = 20)
expression min median
1 factor(y) 2.02ms 2.07ms
2 factor_(y) 460.00µs 485.00µs
simple_mark(base_colsums(df), col_na_counts(df), reps = 20)
expression min median
1 base_colsums(df) 255.00µs 265.00µs
2 col_na_counts(df) 15.00µs 15.00µs
simple_mark(base_seqs(), cheapr_seqs(), reps = 20)
expression min median
1 base_seqs() 1.92ms 1.98ms
2 cheapr_seqs() 30.00µs 40.00µs
simple_mark(
+ do.call(c, list_of_objs),
+ do.call(c_, list_of_objs),
+ c_(.args = list_of_objs),
+ reps = 15
+ )
expression min median
1 do.call(c, list_of_objs) 435.00µs 470.00µs
2 do.call(c_, list_of_objs) 1.67ms 1.75ms
3 c_(.args = list_of_objs) 115.00µs 120.00µs
simple_mark(shallow_copy(df2), deep_copy(df2), reps = 30)
expression min median
1 shallow_copy(df2) 0ns 0ns
2 deep_copy(df2) 0ns 10.00µs結果を表にまとめると次のとおりです。
| 比較対象 | 処理内容 | 実行時間（中央値） | speedup（目安） |
|---|---|---|---|
| cut(x, 10) → cheapr_cut(x, 10) | 10万件の数値を10階級に分割 | 2.34ms → 1.30ms | 約1.8倍 |
| factor(y) → factor_(y) | 10万件の文字ベクトルを因子化 | 2.07ms → 485µs | 約4.3倍 |
| colSums(is.na(df)) → col_na_counts(df) | 100列×1000行の列ごとNA集計 | 265µs → 15µs | 約17.7倍 |
| lapply(seq) → seq_(as_list=TRUE) | 200本の数列を一括生成 | 1.98ms → 40µs | 約49.5倍 |
| do.call(c_, list) | c_()をdo.call経由で呼んだ場合 | 1.75ms | base Rより約3.7倍遅い |
| c_(.args = list) | c_()を.args引数で直接呼んだ場合（推奨） | 120µs | 約3.9倍 |
| shallow_copy(df) → deep_copy(df) | 属性複製の有無によるコピー | 0ns → 10µs | Windowsの時計分解能の影響で0nsとなった可能性あり。参考値 |
＜おすすめのRに関する書籍です＞
Rで楽しむベイズ統計入門［しくみから理解するベイズ推定の基礎］ Data Science Library
ベイズ推定のしくみを、仕組みの理解から丁寧に追う入門書。
www.amazon.co.jp
この記事が誰かの役に立ちますように。
スポンサーリンク
