「統計学関連なんでもあり」の過去ログ--- 043

No.12587　行列の処理では行と列どちらを優先させるべき？　　【波音】　2010/05/07(Fri) 19:09

No.12576に関するスレッドを見ていてふと思ったのですが，行列の扱いで例えば：

f <- function() {
	x <- matrix(numeric(1000*1000), ncol=1000)
	d <- 1
	for(i in 1:1000) {
		for(j in 1:1000) {
			x[i, j] <- d
			d <- d + 1
		}
	}
	return(x)
}

のようなx[1,1] <- 1, x[1,2] <- 2, ...と列方向へ処理を進めていく場合と，これに対して2つ目のfor()のところをx[j, i] <- dとしてx[1,1] <- 1, x[2,1], ...と行方向へ処理を進めていく場合とではどちらが高速なのでしょうか？

system.time(f())を試してみると　x[j, i] <- d　として行方向へと処理を進めていくほうが若干速いようですが・・・

(ちなみにapplyとかを使ったほうがより高速とかいうのはナシで ^_^; )

No.12589　Re: 行列の処理では行と列どちらを優先させるべき？　　【青木繁伸】　2010/05/07(Fri) 21:14

C を初め，多くのプログラム言語では，配列のアクセスは右側の添え字が最も速く動く方向です。ところが，FORTRAN はこれとは逆のアクセスをするのです。まだ世の中に C などなかったころ，FORTRAN が隆盛を誇っており，小さなメモリーと遅いコンピュータでいかに効率的に計算をするかのノウハウがいろいろありました。その一つが，メモリの順にアクセスする（そうしないと，メモリと磁気ディスクの間でデータのスワップが頻発する）というものでした。R の関数などもまだ FORTRAN で書かれたプログラム（実際はそれを C に変換して実行することもある）も多いですね。
R がどちらかはやってみるとよいです。
sim <- function(n)
{
	gc(); gc()
	x <- matrix(0, n, n)
	a <- system.time(for (i in 1:n) for (j in 1:n) x[i, j] <- 1)
	b <- system.time(for (j in 1:n) for (i in 1:n) x[i, j] <- 1)
	return(unname(c(a[1], b[1])))
}
n <- c(100, 500, 1000, 2000, 3000, 4000, 5000)
ans <- sapply(n, sim)

結論は，余り変わらないということかなあ。

行列数サイズ   100   500  1000   2000   3000   4000   5000
行方向       0.036 0.778 3.181 12.538 28.409 51.106 79.759
列方向       0.030 0.793 3.159 12.660 28.550 50.838 79.308
なお，apply，sapply などの速度は for と殆ど変わりません。ベクトル化された計算が速いというのと混同されている気配があります。
> m <- 10000
> x <- numeric(m)
> system.time(for (i in 1:m) x[i] <- sqrt(i))
   ユーザ   システム       経過  
      0.05       0.00       0.05 
> system.time(sapply(1:m, function(i) x[i] <- sqrt(i)))  # これは正しくないし，驚異的に時間を消費する
   ユーザ   システム       経過                          # m が大きくなると，とんでもないことになる
     0.732      0.169      0.889 
> system.time(sapply(1:m, function(i) x[i] <<- sqrt(i))) # 大域変数への付値はお勧めできない
   ユーザ   システム       経過  
     0.068      0.000      0.068 
> system.time(x <- sapply(1:m, sqrt)) # このようにするのではなく
   ユーザ   システム       経過  
      0.02       0.00       0.02 
> system.time(x <- sqrt(1:m))         # このようにすべし
   ユーザ   システム       経過  
     0.000      0.000      0.001      # ほとんど 0

No.12591　Re: 行列の処理では行と列どちらを優先させるべき？　　【MaskedStat】　2010/05/08(Sat) 02:30

SASで以前試してみたことがあり，その頃の記録を引っ張ってみると。。。
500×500行列で，
行方向 1.15 sec
列方向 1.17 sec　でした。
(絶対値が大きいのは当時のPCスペックということでなしで)

C やFORTRANのレベルまでいけば，どのような配列のデータでも解析可能ですが，R/S-PLUSやSASのようなパッケージソフトになると解析用関数によってLONG型，WIDE型のいずれかでないと対応できなくなるので，計算言語の性質が効いてくるように思えます。

No.12593　Re: 行列の処理では行と列どちらを優先させるべき？　　【波音】　2010/05/08(Sat) 09:51

回答ありがとうございます。

う～ん，なるほど。現在でも大規模なデータを扱う場合などは高速化のための書き方を工夫する必要があるのでしょうが，ごく一般的な処理ではコンピュータの性能が格段に良くなっているので，メモリの使われ方など考えなくてもできてしまうということなのだろうと思いました(^_^;)　とはいえ

> system.time(sapply(1:m, function(i) x[i] <- sqrt(i)))

のような恐ろしく時間のかかってしまうやり方もあるわけですから，得られる結果が同じでもどのような書き方が良いのかを意識しておくことは重要なのですね。

>　R/S-PLUSやSASのようなパッケージソフトになると解析用関数によってLONG型，WIDE型のいずれかでないと対応できなくなるので，計算言語の性質が効いてくるように思えます。

奥が深いことですね。やはりFORTRANの時代にプログラミングを学んだ人たちは，おそらくコンピュータの仕組みを理解していないと難しかったのでしょうが，現在は極端な話，仕組みなど知らずとも，コードの書式さえ覚えてしまえばある程度はこなせてしまうのだから(私もその1人ですが)便利というかなんというか。。。

No.12594　Re: 行列の処理では行と列どちらを優先させるべき？　　【波音】　2010/05/08(Sat) 10:15

度々すいません(^_^;)

ベクトル単位で処理を行うのが好ましいというのは，(最も単純な例で)例えば：

> x <- numeric(100000)
> for(i in 1:100000) x[i] <- i

というコードだと繰り返しの度に(i回)オブジェクトxにアクセスするが，

> y <- numeric(100000)
> y <- 1:100000

ならば1回オブジェクトyにアクセスするだけで済むから，こちらの方が"好ましい"ということでよいのでしょうか？

No.12595　Re: 行列の処理では行と列どちらを優先させるべき？　　【青木繁伸】　2010/05/08(Sat) 20:18

その例も当てはまるでしょう。なお，後者は，y <- numeric(100000) は不要で，y <- 1:100000 だけでよいですね。

No.12597　Re: 行列の処理では行と列どちらを優先させるべき？　　【波音】　2010/05/08(Sat) 22:39

ありがとうございました＾＾