【数値情報処理a】第2回 数の表現と誤差
前回の課題
2進数と16進数
コンピュータ内部では、データ（数値、文字、画像など）は全て2 進数の状態で保存している。
2 進数の一桁の情報量をビット（bit）という。 さらに、8 ビットの情報量をバイト（byte）という。
人間にとっては普段使っている 10 進数の方が 2 進数より扱いやすいが、2 進数と 10 進数では桁の対応がシンプルではない。 そこで2 進数四桁（4 ビット）を一桁で表すことができる16 進数がよく使われる。
16 進数では一つの桁を $0$ から $9$ までの 10 文字だけでは表せない。 したがって a から f までのアルファベット 6 文字を加えて合わせて16文字で一つの桁を表す。
| 2進数 | 10進数 | 16進数 |
|---|---|---|
| 0 | 0 | 0 |
| 1 | 1 | 1 |
| 10 | 2 | 2 |
| 11 | 3 | 3 |
| 100 | 4 | 4 |
| 101 | 5 | 5 |
| 110 | 6 | 6 |
| 111 | 7 | 7 |
| 1000 | 8 | 8 |
| 1001 | 9 | 9 |
| 1010 | 10 | a |
| 1011 | 11 | b |
| 1100 | 12 | c |
| 1101 | 13 | d |
| 1110 | 14 | e |
| 1111 | 15 | f |
| 10000 | 16 | 10 |
| … | … | … |
例えば 10 進数の整数で $51966$ は、2 進数と 16 進数では以下のように表される。 \[ 51966=(1100\ 1010\ 1111\ 1110)_{2}=(\mathrm{cafe})_{16} \] かっこの下付きの数字は、かっこの中の数字の進数を表す。
N進数の整数
N進数の整数から10進数への変換
\(m\)桁の N 進数の整数 $(a_{m}a_{m-1}\cdots a_{3}a_{2}a_{1})_{N}$ から 10 進数への変換は以下で計算する。
\begin{eqnarray} &&\sum_{i=1}^{m}a_{i}\times N^{i-1}\\ &=&a_{m}\times N^{m-1}+a_{m-1}\times N^{m-2}+\cdots+a_{3}\times N^{2}+a_{2}\times N^{1}+a_{1}\times N^{0} \end{eqnarray}
変換例1
16 進数の整数 $(\mathrm{cafe})_{16}$ を 10 進数に変換する。
まず、各桁を 10 進数に変換すると以下のようになる。 \begin{eqnarray} a_{4}&=&(\mathrm{c})_{16}=12\\ a_{3}&=&(\mathrm{a})_{16}=10\\ a_{2}&=&(\mathrm{f})_{16}=15\\ a_{1}&=&(\mathrm{e})_{16}=14 \end{eqnarray} これより \begin{eqnarray} &&12\times 16^{3}+10\times 16^{2}+15\times 16^{1}+14\times 16^{0}\\ &=&12\times 4096+10\times 256+15\times 16+14\\ &=&51966 \end{eqnarray} となる。
変換例2
2 進数の整数 $(1100\ 1010\ 1111\ 1110)_{2}$ を 10 進数に変換する。 \begin{eqnarray} &&1\times 2^{15}+1\times 2^{14}+0\times 2^{13}+0\times 2^{12}+\cdots+1\times 2^{3}+1\times 2^{2}+1\times 2^{1}+0\times 2^{0}\\ &=&32768+16384+2048+512+128+64+32+16+8+4+2\\ &=&51966 \end{eqnarray} となる。
10進数の整数からN進数への変換
10 進数の整数を N 進数に変換するには $N$ で割った余りを使う。 $N$ で割って商と余りを求め、その商をさらに $N$ で割って商と余りを求めることを繰り返す。 商が $0$ になると終了する。
この原理を簡単に説明するために整数 $x$ が次の形に表せたとしよう。 \[ x = a_{3}\times N^{2}+a_{2}\times N^{1}+a_{1}\times N^{0} \] ここで $a_{3},a_{2},a_{1}$ はそれぞれ $0$ から $N-1$ の間の整数である。
$N$ で割って商と余りを求め、その商をさらに $N$ で割って商と余りを求めることを繰り返した結果は以下である。
| 商 | 余り | ||
|---|---|---|---|
| $(a_{3}\times N^{2}+a_{2}\times N^{1}+a_{1}\times N^{0})\div N=$ | $a_{3}\times N^{1}+a_{2}\times N^{0}$ | $\cdots$ | $a_{1}$ |
| $(a_{3}\times N^{1}+a_{2}\times N^{0})\div N=$ | $a_{3}\times N^{0}$ | $\cdots$ | $a_{2}$ |
| $(a_{3}\times N^{0})\div N=$ | $0$ | $\cdots$ | $a_{3}$ |
最後は商が $0$ になったので終了した。
余りを下から上に並べたものが整数 $x$ の N 進数である。 \[ x=(a_{3}a_{2}a_{1})_{N} \]
変換例
10進数の $6$ を 2 進数に変換する。
| 商 | 余り | ||
|---|---|---|---|
| $6\div 2=$ | $3$ | $\cdots$ | $0$ |
| $3\div 2=$ | $1$ | $\cdots$ | $1$ |
| $1\div 2=$ | $0$ | $\cdots$ | $1$ |
余りを下から上に並べると $6$ の 2 進数が得られる。 \[ 6=(110)_{2} \] 確認のために、この 2 進数を 10 進数に変換してみると、確かに 6 になることが分かる。 \[ 1\times 2^{2}+1\times 2^{1}+0\times 2^{0}=4+2+0=6 \]
RにおけるN進数の整数
N進数の整数から10進数への変換
N 進数の整数から10進数への変換は上述の式で計算すればよい。
例えば、16 進数の整数 $(\mathrm{cafe})_{16}$ から 10 進数への変換は次のようになる。
> 12 * 16 ^ 3 + 10 * 16 ^ 2 + 15 * 16 ^ 1 + 14 * 16 ^ 0 [1] 51966
R ではもっと簡単な方法として strtoi(文字列, N) で変換できる。
> strtoi("cafe", 16) [1] 51966 > strtoi("1100101011111110", 2) [1] 51966
N 進数の数は文字列として入力する必要がある。
R では「"」（ダブルクォーテーションマーク）で囲ったものは文字列として扱われる。
10進数の整数から16進数への変換
10 進数から 16 進数への変換は sprintf() か as.hexmode() を使う。
> sprintf("%x", 51966) [1] "cafe" > as.hexmode(51966) [1] "cafe"
10 進数から 2 進数への変換は intToBits() があるが、使いづらいので省略する。
小テスト①
Moodle Server（非公式）で第2回の小テスト①を受験しなさい。
コンピュータにおける整数の表現
コンピュータで整数は符号部（$\pm$：プラスかマイナス）と正の整数で表現される。
\[ 整数=\pm(正の整数) \]
整数は無限に存在するので、有限な容量しか持たないコンピュータでは全ての整数を表現することはできない。 そこで、整数にある決まったデータ長を割り当てて、十分大きな有限の範囲の整数が扱えるようにしてある。
R では一つの整数に32ビット（4バイト）のデータ長が使われる。 各ビットは以下のように割り当てられる。
- 符号部：1 ビット
- 正の整数：31 ビット
各ビットの格納のイメージは下図である1)。
ここで $a_{i}$（$i=1,2,\cdots,31$）はそれぞれ 0 か 1 のどちらかである（2進数）。
\[ 正の整数=\sum_{i=1}^{31}a_{i}\times 2^{i-1}=a_{1}\times 2^{0}+a_{2}\times 2^{1}+\cdots+a_{31}\times 2^{30} \]
正の整数の最大値
正の整数の最大値は $a_{i}$ がすべて $1$ の場合なので \[ 正の整数の最大値=\sum_{i=1}^{31}1\times 2^{i-1}=\frac{1-2^{31}}{1-2}=2^{31}-1=2147483647 \] になる。 ここで、等比数列（初項 $1$、公比 $2$）の和の公式を使った。
正の整数の最小値
正の整数の最小値は $a_{i}$ がすべて $0$ の場合なので \[ 正の整数の最小値=\sum_{i=1}^{31}0\times 2^{i-1}=0 \] になる。
整数の範囲
正の整数の最大値より、R で表現可能な整数の範囲は \[ -2147483647\le 32ビット整数\le 2147483647 \] になる。
Rにおける整数
R では数字は一般的に実数型 “numeric” として扱われるため、たとえ整数であっても実数とみなされる。
強制的に整数型 “integer” として扱うには数字の最後に L を付けるか as.interger(数字) を使う必要がある。
> class(123) [1] "numeric" > class(123L) [1] "integer" > class(as.integer(123)) [1] "integer"
整数型 “integer” での最大値は $2147483647$ である。 この値を超えると自動的に “numeric” になる。
> class(2147483647L) [1] "integer" > class(2147483648L) [1] "numeric" 警告メッセージ: 非整数値 2147483648L に L がついています; 実数値を使用します
N進数の小数
N進数の小数から10進数への変換
\(m\)桁の N 進数の小数 $(0.b_{1}b_{2}b_{3}\cdots b_{m-1}b_{m})_{N}$ から 10 進数への変換は以下で計算する。
\begin{eqnarray} &&\sum_{i=1}^{m}b_{i}\times N^{-i}\\ &=&b_{1}\times N^{-1}+b_{2}\times N^{-2}+b_{3}\times N^{-3}+\cdots+b_{m-1}\times N^{-m+1}+b_{m}\times N^{-m} \end{eqnarray}
変換例
2 進数の小数 $(0.101)_{2}$ から 10 進数に変換する。 \[ 1\times 2^{-1}+0\times 2^{-2}+1\times 2^{-3}=\frac{1}{2}+\frac{1}{8}=0.5+0.125=0.625 \] となる。
10進数の小数から2進数への変換
10 進数で $1$ より小さい小数を 2 進数に変換するには $2$ をかけた整数部を使う。 $2$ をかけて整数部と小数部に分け、その小数部にさらに $2$ をかけて整数部と小数部に分けることを繰り返す。 小数部が $0.0$ になると終了する。
この原理を簡単に説明するために小数 $y$ （$0\le y < 1$）が次の形に表せたとしよう。 \[ y = b_{1}\times 2^{-1}+b_{2}\times 2^{-2}+b_{3}\times 2^{-3} \] ここで、$b_{1},b_{2},b_{3}$ はそれぞれ $0$ か $1$ のどちらかである。
$2$ をかけて整数部と小数部に分け、その小数部にさらに $2$ をかけて整数部と小数部に分けることを繰り返した結果は以下である。
| 整数部 | 小数部 | ||
|---|---|---|---|
| $(b_{1}\times 2^{-1}+b_{2}\times 2^{-2}+b_{3}\times 2^{-3})\times 2=$ | $b_{1}$ | $+$ | $b_{2}\times 2^{-1}+b_{3}\times 2^{-2}$ |
| $(b_{2}\times 2^{-1}+b_{3}\times 2^{-2})\times 2=$ | $b_{2}$ | $+$ | $b_{3}\times 2^{-1}$ |
| $(b_{3}\times 2^{-1})\times 2=$ | $b_{3}$ | $+$ | $0.0$ |
最後は小数部が $0.0$ になったため終了した。
この整数部を上から下に並べたものが小数 $y$ の 2 進数である。 \[ y = (0.b_{1}b_{2}b_{3})_{2} \]
変換例1
10 進数の小数 $0.625$ を 2 進数に変換する。
| 整数部 | 小数部 | ||
|---|---|---|---|
| $0.625\times 2=$ | $1$ | $+$ | $0.25$ |
| $0.25\times 2=$ | $0$ | $+$ | $0.5$ |
| $0.5\times 2=$ | $1$ | $+$ | $0.0$ |
整数部を上から下に並べると $0.625$ の 2 進数が得られる。 \[ 0.625=(0.101)_{2} \] この 2 進数を 10 進数に変換してみると、確かに $0.625$ になることが分かる。 \[ 1\times 2^{-1}+0\times 2^{-2}+1\times 2^{-3}=0.625 \]
変換例2
10 進数の小数 $0.3$ を 2 進数に変換する。
| 整数部 | 小数部 | ||
|---|---|---|---|
| $0.3\times 2=$ | $0$ | $+$ | $0.6$ |
| $0.6\times 2=$ | $1$ | $+$ | $0.2$ |
| $0.2\times 2=$ | $0$ | $+$ | $0.4$ |
| $0.4\times 2=$ | $0$ | $+$ | $0.8$ |
| $0.8\times 2=$ | $1$ | $+$ | $0.6$ |
| $0.6\times 2=$ | $1$ | $+$ | $0.2$ |
| $0.2\times 2=$ | $0$ | $+$ | $0.4$ |
| $0.4\times 2=$ | $0$ | $+$ | $0.8$ |
| $0.8\times 2=$ | $1$ | $+$ | $0.6$ |
| $\vdots$ | $\vdots$ | $+$ | $\vdots$ |
すぐに分かるが、この例では手順をいくら繰り返しても小数部が $0.0$ になることはないため、終わりがない。
整数部を上から下に並べると $0.3$ の 2 進数が得られる。 途中から $1001$ が無限に繰り返し現れる。
\[ 0.3=(0.0\color{blue}{1001}1001\color{blue}{1001}1001\color{blue}{1001}1001\color{blue}{1001}1001\color{blue}{1001}1001\color{blue}{1001}1001\color{blue}{1001}\cdots)_{2} \]
10 進数であっても $\frac{1}{3}=0.33333\cdots$ のように無限に同じ数を繰り返す小数（循環小数）があるので、N 進数でそのような小数があってもおかしくはない。
N進数の桁の上げ下げ
N 進数の数の桁を 1 つ上げるには $N$ をかける。 \[ (a_{4}a_{3}a_{2}a_{1})_{N}\times N=(a_{4}a_{3}a_{2}a_{1}0)_{N} \] N 進数の数の桁を 2 つ上げるには $N$ を 2 回かける（$N^{2}$ をかける）。 \[ (a_{4}a_{3}a_{2}a_{1})_{N}\times N^{2}=(a_{4}a_{3}a_{2}a_{1}00)_{N} \]
逆に N 進数の数の桁を 1 つ下げるには $N$ で割る。 \[ (a_{4}a_{3}a_{2}a_{1})_{N}\div N=(a_{4}a_{3}a_{2}.a_{1})_{N} \] N 進数の数の桁を 2 つ下げるには $N$ で 2 回割る（$N^{2}$ で割る）。 \[ (a_{4}a_{3}a_{2}a_{1})_{N}\div N^{2}=(a_{4}a_{3}.a_{2}a_{1})_{N} \]
10進数の例
| 10進数 |
|---|
| $1234\times 10=12340$ |
| $1234\times 10^{2}=123400$ |
| $1234\div 10=123.4$ |
| $1234\div 10^{2}=12.34$ |
2進数の例
| 2進数 | 10進数との対応 |
|---|---|
| $(1001)_{2}\times 2=(10010)_{2}$ | $9\times 2=18$ |
| $(1001)_{2}\times 2^{2}=(100100)_{2}$ | $9\times 2^{2}=36$ |
| $(1001)_{2}\div 2=(100.1)_{2}$ | $9\div 2=4.5$ |
| $(1001)_{2}\div 2^{2}=(10.01)_{2}$ | $9\div 2^{2}=2.25$ |
小テスト②
Moodle Server（非公式）で第2回の小テスト②を受験しなさい。
コンピュータにおける実数の表現
実数も無限に存在するので、コンピュータでは全ての実数を表現することはできない。 実数でもある決まったデータ長を割り当てて、有限な範囲の実数が扱えるようにしてある。
コンピュータ内部では実数も 2 進数で保存される。 実数は浮動小数点数で実現され、符号部（$\pm$）、指数部、仮数部を使って表現される。 一般的に使われている IEEE 7542)では以下のように表現される。 \[ 浮動小数点数=\pm(1+仮数部)\times 2^{(指数部-E_{0})} \] ここで $0\le 仮数部 < 1$ である。 $2$ の指数にある $E_{0}$ は定数で、その値は精度に応じて定められている。
浮動小数点数の精度は、使われるデータ長によって単精度、倍精度、四倍精度などといわれる。 四倍長精度以上を多倍長精度ということがある。
| 精度 | データ長 | $E_{0}$ |
|---|---|---|
| 単精度 | 32 ビット | $127$ |
| 倍精度 | 64 ビット | $1023$ |
| 四倍精度 | 128 ビット | $16383$ |
R では実数に倍精度浮動小数点数が使われる。 各ビットは以下のように割り当てられる。
- 符号部：1 ビット
- 指数部：11 ビット
- 仮数部：52 ビット
各ビットの格納のイメージは下図である3)。
ここで $b_{i}$ はそれぞれ $0$ か $1$ のどちらかである（2進数）。
仮数部は 10 進数で以下のように定義される。 \[ 仮数部=\sum_{i=1}^{52}b_{i}\times 2^{-i} \]
仮数部の最大値
仮数部の最大値は $b_{i}$ が全て $1$ の場合なので \[ 仮数部の最大値=\sum_{i=1}^{52}2^{-i}=\sum_{i=1}^{52}2^{-1}\times (2^{-1})^{i-1}=\frac{2^{-1}(1-2^{-52})}{1-2^{-1}}=1-2^{-52} \] になる。 ここで、等比数列（初項 $2^{-1}$、公比 $2^{-1}$）の和の公式を使った。
仮数部の最小値
仮数部の最小値は $b_{i}$ が全て $0$ の場合なので \[ 仮数部の最小値=0 \] になる。
仮数部の範囲
仮数部の最大値と最小値より、仮数部の範囲は \[ 0\le 仮数部\le 1-2^{-52} \] になる。
指数部の範囲
指数部は 11 ビットの正の整数なので $0$ から $2047$ の範囲になるが、$0$ には特殊な意味を持たせるので（後述）、 利用できる範囲は \[ 1\le 指数部\le 2047 \] である。 倍精度では $E_{0}=1023$ が使われるので、$2$ の指数の範囲は \[ -1022\le 指数部-E_{0}\le 1023 \] になる。
正の浮動小数点数の範囲
正の浮動小数点数で表現可能な範囲は $0$ を除いて \[ (1+0)\times 2^{-1022}\le 正の浮動小数点数 \le (2-2^{-52})\times 2^{1023} \] になり、計算して四捨五入すると \[ 2.225074\times 10^{-308}\le 正の浮動小数点数 \le 1.797693\times 10^{308} \] になる。
ただし、この範囲では実数が $0$ の場合を表現できない。 そこで、仮数部と指数部がともに $0$ の場合を実数の $0$ と定義する。
浮動小数点数の例
例として 10 進数の小数 $9.3$ を浮動小数点数にしてみよう。
実数の2進数への変換
まず、$9.3$ を 2 進数に変換する。
① 整数部と小数部に分離
$9.3$ の整数部は $9$、小数部は $0.3$ である。
② 整数部を 2 進数へ変換
整数部 $9$ を 2 進数に変換する。 \[ 9=(1001)_{2} \]
③ 小数部を2進数へ変換
小数部 $0.3$ を 2進数に変換する。
\[ 0.3=(0.0\color{blue}{1001}1001\color{blue}{1001}1001\color{blue}{1001}1001\color{blue}{1001}1001\color{blue}{1001}1001\color{blue}{1001}1001\color{blue}{1001}\cdots)_{2} \]
以上より、$9.3$ は 2 進数で次のように表せる。 \[ 9.3=(1001.0\color{blue}{1001}1001\color{blue}{1001}1001\color{blue}{1001}1001\color{blue}{1001}1001\color{blue}{1001}1001\color{blue}{1001}1001\color{blue}{1001}\cdots)_{2} \] これは無限に続く数なので、どこかで丸めてコンピュータで扱えるようにしたのが浮動小数点数である。
浮動小数点数への変換
次にこれを浮動小数点数にする。
まず、2 進数の小数の桁を 3 桁ずらして $(1+仮数部)$ になるようにする。これを正規化という。 \begin{align*} 9.3&=(1001.0\color{blue}{1001}1001\color{blue}{1001}1001\color{blue}{1001}1001\color{blue}{1001}1001\color{blue}{1001}1001\color{blue}{1001}1001\color{blue}{1001}\cdots)_{2}\\ &=(1.0010\color{blue}{1001}1001\color{blue}{1001}1001\color{blue}{1001}1001\color{blue}{1001}1001\color{blue}{1001}1001\color{blue}{1001}1001\color{blue}{1001}\cdots)_{2}\times 2^{3}\\ &=(1+0.0010\color{blue}{1001}1001\color{blue}{1001}1001\color{blue}{1001}1001\color{blue}{1001}1001\color{blue}{1001}1001\color{blue}{1001}1001\color{blue}{1001}\cdots)_{2}\times 2^{3} \end{align*} これを \[ 浮動小数点数=\pm(1+仮数部)\times 2^{(指数部-E_{0})} \] と比べると \begin{align*} 指数部&=3+E_{0}=1026\\ 仮数部&=(0.0010\color{blue}{1001}1001\color{blue}{1001}1001\color{blue}{1001}1001\color{blue}{1001}1001\color{blue}{1001}1001\color{blue}{1001}1001\color{blue}{1001}\cdots)_{2} \end{align*} である。 仮数部は倍精度で52ビットまでの数字しか扱えないので、 仮数部の小数点第53位を零捨一入4)で丸めて小数点第52位までにする。 \begin{align*} 仮数部&=(0.0010\color{blue}{1001}1001\color{blue}{1001}1001\color{blue}{1001}1001\color{blue}{1001}1001\color{blue}{1001}1001\color{blue}{1001}1001\color{blue}{1}\!\!\!\!\!\!\!\!\!\backslash\backslash\backslash\cdots)_{2}\\ &\simeq(0.0010\color{blue}{1001}1001\color{blue}{1001}1001\color{blue}{1001}1001\color{blue}{1001}1001\color{blue}{1001}1001\color{blue}{1001}1010)_{2} \end{align*}
したがって、9.3 は 2 進数の浮動小数点数で \begin{align*} 9.3&\simeq(1 + 0.0010\color{blue}{1001}1001\color{blue}{1001}1001\color{blue}{1001}1001\color{blue}{1001}1001\color{blue}{1001}1001\color{blue}{1001}1010)_{2}\times 2^{3}\\ &\simeq(1.0010\color{blue}{1001}1001\color{blue}{1001}1001\color{blue}{1001}1001\color{blue}{1001}1001\color{blue}{1001}1001\color{blue}{1001}1010)_{2}\times 2^{3}\\ &\simeq(1001.0\color{blue}{1001}1001\color{blue}{1001}1001\color{blue}{1001}1001\color{blue}{1001}1001\color{blue}{1001}1001\color{blue}{1001}1010)_{2} \end{align*} となる。
Rにおける浮動小数点数
R で直接 $9.3$ を入力してみよう。
> 9.3 [1] 9.3
ちゃんと $9.3$ が表現されているように見えるが、これは R ではデフォルトで 7 桁までの表示になっているためである。 表示桁数を最大の 22 桁に変更してみると
> print(9.3, digits = 22) # 表示桁数を22桁に変更 [1] 9.300000000000000710543
と表示され、浮動小数点数では $9.3$ を正確に表現できていないことが分かる。 $9.3$ をこちらが入力しても、コンピュータ内部では上述したように浮動小数点数による丸めが行われているので、$9.3$ ぴったりにはならない。
この丸めによって生じる誤差を丸め誤差という。
先ほど求めた $9.3$ の 2 進数の浮動小数点数 \[ (1001.0\color{blue}{1001}1001\color{blue}{1001}1001\color{blue}{1001}1001\color{blue}{1001}1001\color{blue}{1001}1001\color{blue}{1001}1010)_{2} \] を 10 進数に変換して丸め誤差を確かめよう。
次のプログラム bin2dec.R は 2 進数の実数を 10 進数に変換するプログラムである。
x に 2 進数の浮動小数点数を文字列で入力する必要がある。
Rmpfr は任意精度用のパッケージで、70 桁の精度まで丸めずに計算するために使っている。
- bin2dec.R
- # 2進数から10進数への変換
- # 任意精度用のパッケージ
- # 精度（ビット）
- pbit <- 70
- # 2進数（文字列）
- x <- "1001.0100110011001100110011001100110011001100110011010"
- # xを整数部と小数部に分ける
- # 整数部の桁数
- # 整数部と小数部を結合して1文字ずつに分割
- # 1桁ずつ計算して足し合わせ
- r <- 0
- r <- mpfr(r + strtoi(a[i], 2) * 2 ^ (idigit - i), pbit)
- }
- # 10進数の出力
結果
[1] 9.300000000000000710543
「print(9.3, digits = 22)」と同じ結果が得られることが分かる。
正の浮動小数点数の最大値
R で扱うことができる正の浮動小数点数の最大値、つまり R でどこまで大きい数を表現できるかは .Machine$double.xmax で確認できる。
> .Machine$double.xmax [1] 1.797693e+308
1.797693e+308 は $1.797693\times 10^{308}$ のことである。
1.797693e+308 を超えると Inf（無限大）として扱われる。
> 1.797693e+308 [1] 1.797693e+308 > 1.797694e+308 [1] Inf
正の浮動小数点数の最小値
R で扱うことができる正の浮動小数点数の（$0$ より大きい）最小値、つまり R でどこまで数を小さく分けられるかは .Machine$double.xmin で確認できる。
> .Machine$double.xmin [1] 2.225074e-308
2.225074e-308 は $2.225074\times 10^{-308}$ のことである。
数値計算で生じる誤差
- 丸め誤差
- 打ち切り誤差
- 桁落ち
- 情報落ち
丸め誤差
仮数部の上限の桁数を超えたときに丸めることによって生じる誤差である。
打ち切り誤差
無限に繰り返し行って結果が得られる計算処理を途中で打ち切ることにより生じる誤差である。
例えば、等比数列の和 \[ S_{n}=\sum_{i=1}^{n}(0.5)^{i-1} \] の $n\rightarrow\infty$ の極限を考えよう。
等比数列の和の公式から \[ \lim_{n\rightarrow\infty}S_{n}=\lim_{n\rightarrow\infty}\frac{1-(0.5)^{n}}{1-0.5}=\frac{1}{0.5}=2 \] となることは分かる。
R で $S_{n}$ を計算したらどうなるか見てみよう。 無限回の和をとることは不可能なので $n=10$ で近似して計算すると次のようになる。
> s <- 0; for (i in 1 : 10) {s <- s + 0.5 ^ (i - 1)}; print(s) [1] 1.998047
和を 10 回で打ち切っているため、$2$ にならない。 この場合の打ち切り誤差は $2-1.998047=0.001953$ である。
$n$ を増やせば打ち切り誤差は小さくなっていくことが確認できる。
桁落ち
値がほぼ等しく丸め誤差を持つ二つの数値の引き算を行った時に生じる誤差のことである。
有名な例は二次方程式の解の計算で起こる桁落ちである。
$ax^{2}+bx+c=0$ の解は、解の公式より \[ x_{\pm}=\frac{-b\pm\sqrt{b^{2}-4ac}}{2a} \] であるが、$b$ に比べて $a,c$ が非常に小さいときに片方の解が桁落ちを起こしやすい。
$a=b=1, c=10^{-8}$ として片方の解 $x_{+}$ を計算すると、正しくは \[ x_{+}= -1.00000001000\cdots\times 10^{-8} \] となるが、実際に R で計算すると
> a <- 1; b <- 1;c <- 1e-8; print((-b + sqrt(b ^ 2 - 4 * a * c))/(2 * a), digits = 22) [1] -1.000000010575874398455e-08
となり、小数点第10位に誤差が生じる。
この桁落ちを防ぐには、工夫してほぼ等しい数同士の引き算を無くせばよい。 $x_{+}$ は次のように変形できる。 \[ x_{+}=\frac{(-b+\sqrt{b^{2}-4ac})(-b-\sqrt{b^{2}-4ac})}{2a(-b-\sqrt{b^{2}-4ac})}=-\frac{2c}{b+\sqrt{b^{2}-4ac}} \] 変形した式だと、ほぼ等しい数同士の引き算がないので桁落ちを防ぐことができる。
この式で計算すると
> a <- 1; b <- 1;c <- 1e-8; print(-2 * c / (b + sqrt(b ^ 2 - 4 * a * c)), digits = 22) [1] -1.000000010000000354983e-08
となる。
情報落ち
絶対値が大きい数と小さい数を足したり引いたりしたときに小さい数の情報が落ちてしまうことである。
例えば、$10^{16}$ に $10$ を足してみよう。
> print(1e16 + 10, digits = 22) [1] 10000000000000010
$10000000000000010$ という結果になり、情報が残っているように見える。
一方、$10^{16}$ に $1$ を足すと情報が落ちてしまうので $10000000000000001$ にならない。
> print(1e16 + 1, digits = 22) [1] 1e+16
$10^{16}$ に $1$ を 10 回足したとしても情報が落ちているので $10000000000000010$ にならない。
> print(1e16 + 1 + 1 + 1 + 1 + 1 + 1 + 1 + 1 + 1 + 1, digits = 22) [1] 1e+16
課題
Moodle Server（非公式） で第2回の課題を行いなさい。
締め切り：2025年4月23日（水）20時