四捨五入関数は round ではない（かも）

三省堂 大辞林によれば，四捨五入とは，「必要とする位の次の位の数字が四以下ならばこれを切り捨て、五以上ならば切り上げて、上の位に一を加える方法」とある。

R では，round 関数がこれを行う。

1.5 を，小数第1位で四捨五入を行うと 2 である。
> round(1.5)
[1] 2

しかし 0.5 を，小数第1位で四捨五入を行うと 0 である。

 > round(0.5)
[1] 0

よって，round は四捨五入を行う関数ではない。

round は，偶数への丸め（最近接丸め，JIS 丸め，ISO 丸め，銀行丸めなどとも呼ばれる）を行う関数である。これは，「端数が 0.5 より小さいなら切り捨て、端数が 0.5 より大きいならは切り上げ、端数がちょうど 0.5 なら切り捨てと切り上げのうち結果が偶数となる方へ丸める。（JIS Z 8401 の規則A）」として定められており、規則B（四捨五入）より「望ましい」とされている。

インターネット上の多くのページでは，プログラム言語（処理系）によって違いがあるとかの記述はあるが，実際に四捨五入を行うにはどのような関数を定義すればよいかはほとんど書かれていない。

そこで，R においては以下の関数を定義する。

round2 = function(x, d=0) {
    p = 10^d
    return((x * p * 2 + 1) %/% 2 / p)
}

この関数はまず，四捨五入を行う位置 d により，対象とする数 x に p = 10^d をかける。たとえば x = 0.45, d = 1 ならば x = 0.45 * 10^1 = 4.5 となる。ここで小数点以下を四捨五入すれば 5 となる。同じことをするのに x = (x * 2 + 1) %/% 2 としている。x = (4.5 * 2 + 1) = 10 を 2 で割って整数部を取る。つまり x = 5 となるので四捨五入していることになる。なお，(x * p * 2 + 1) %/% 2  は as.integer(x * p + 0.5) と同じである。x * p の整数部 と 0.5 は共にコンピュータ上で正確に表現できる数なので，誤差はない。最後に得られた数を p で割る。つまり x = 5 / 10 = 5 が答えである。

> round2(0.45, 1)
[1] 0.5

これに対して，元からある round は以下のような結果になる。
> round(0.45, 1)
[1] 0.4

解釈は，「0.45 という数値がコンピュータ上で正確に表されているなら，小数点第2位を丸めるときに，切り上げて 0.5 とする場合も，切り下げて 0.4 とする場合もその差の絶対値は 0.05 である。このように丸めの候補値が 2 つあるときは，結果の末尾が偶数のほうを採用する。つまり，この場合だと 0.4 になる」ということである。

注意すべきは「0.45 という数値がコンピュータ上で正確に表されているなら」というところ。0.45 は 10 進数では循環小数ではないが 2 進数では循環小数で内部的には誤差を含む数値であることもある。

0.45 ≒ 0x1.ccccccccccccdp-2 = (1+12/16^1 + 12/16^2 + … + 12/16^12 + 13/16^13) *2 ^-2

> round(0.8500000000000000, 1)
[1] 0.8

となるのは正しい結果であるが，

> round(0.8500000000000001, 1)
[1] 0.8

は 0.8500000000000001 が内部的に正しく表現されているなら，本来は 0.9 にならないといけない。

> round(0.8500000000000002, 1)
[1] 0.9

は正しい結果である。

つまり，コンピュータ内部では，0.8500000000000001 と 0.8500000000000000 は異なる数であるが，

> 0.8500000000000001 == 0.8500000000000000
[1] FALSE

その差は末尾 1 ビット

> 0.8500000000000001 - 0.8500000000000000 < EPSILON
[1] TRUE

> sprintf("%a %a %a %a %a", EPSILON, 0.8499999999999999, 0.8500000000000000, 0.8500000000000001, 0.8500000000000002)

EPSILON    =　0x1p-52
0.8499999999999999　=　0x1.b333333333332p-1
0.8500000000000000　=　0x1.b333333333333p-1
0.8500000000000001　=　0x1.b333333333334p-1
0.8500000000000002　=　0x1.b333333333335p-1"

いずれの場合も，四捨五入関数 round2 は正しく機能している。

> round(0.8499999999999999, 1); round2(0.8499999999999999, 1)
[1] 0.8
[1] 0.8
> round(0.85, 1); round2(0.85, 1)
[1] 0.8
[1] 0.9
> round(0.8500000000000001, 1); round2(0.8500000000000001, 1)
[1] 0.8
[1] 0.9
> round(0.8500000000000002, 1); round2(0.8500000000000002, 1)
[1] 0.9
[1] 0.9