SAS 數值儲存方式和精度問題

Snoopy1866發表於2022-03-18

本文連結：https://www.cnblogs.com/snoopy1866/p/16021137.html

1 數值儲存方式

SAS使用8個位元組儲存數值，使用浮點計數法表示數值。
浮點計數法由4個部分組成：符號位（Sign）、基數（Base）、指數（Exponent）、尾數（Mantissa），這4個部分分別提供以下作用：

  • 符號位：代表該數值是正數還是負數
  • 基數：SAS系統中的基數預設為2。
  • 指數：代表基數被乘的倍數
  • 尾數：定義數值範圍的一個小數

將計算機內部儲存的數值轉化為10進位制數值的公式為：\(\color{red}{Sign}\)*(\(\color{purple}{Mantissa}\) * Base$\color{green}{Exponent}$)
下圖展示了SAS用8個位元組表示浮點數的具體情況：

即第1位為符號位，第2-12位為指數位，13-64位為尾數位，符號位中，0表示正數，1表示負數。
例如：Sign = 1, Exponent = 3, Mantissa = 1492，在10進位制下代表：(-1) * (0.1492 * 103) = -149.2，在2進位制下表示：(-1) * (0.1492 * 22) = -1.1936;

2 產生的問題

不是所有浮點數都能用8個位元組精確地表示，某些浮點數甚至無法用有限個位元組精確地表示，這在任何進位制下都是存在的問題。
例如：10進位制下的0.1在16進位制下表示為3FB999999999999999999999999999999...，2/3在10進位制下表示為0.666666666666666666666...
當計算機面對這種情況時，有兩種選擇：
（1）Truncation：2/3+2/3+2/3 = 0.666666 + 0.666666 + 0.666666 = 1.999998
（2）Rounding: 2/3+2/3+2/3 = 0.666667 + 0.666667 + 0.666667 = 2.000001
但這兩種選擇都無法完美解決精度問題，而且隨著運算次數的累計，精度問題會愈發凸顯。

來看下面一個例子：

/*DO 迭代生成0.1~0.9*/
data s1;
    do a = 0.1 to 0.9 by 0.1;
        index + 1;
        output;
    end;
run;

/*手動輸入0.1~0.9*/
data s2;
    do b = 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9;
        index + 1;
        output;
    end;
run;

/*合併s1-s2*/
data s3;
    merge s1 s2;
    by index;
    if a = b then A_eq_B = "Y";
run;

proc print data = s3;
    var index a b A_eq_B;
run;


可以看到PROC PRINT列印出來的結果中顯示，A，B兩列顯示完全一致，但其中0.3,0.8和0.9卻出現了A≠B的情況。
為了能夠看到SAS內部儲存的實際數值，使用十六進位制輸出格式輸出結果：

proc print data = s3;
    var index a b A_eq_B;
    format a hex16. b hex16.;
run;


由此可見，是DO迴圈迭代計算導致的誤差累計，最終導致A和B兩列實際儲存的數值存在微小差異。
這種精度損失可能會在下述場景中導致問題：

  • 使用 IF，SELECT，WHERE語句比較兩個不同計算方法得出的數值，或與顯式指定的數值比較時；
  • 對資料集取子集時；
  • 使用 PROC COMPARE 過程比較兩個資料集時
  • 對基於計算產生的浮點數進行分類分析時
  • PROC REPORT或其他過程步的輸出可能會顯示出奇怪的小數（例如：-0.00）

3 解決辦法

以下3種方式可以解決SAS中存在的大多數精度問題：

  • 使用 ROUND 函式
  • 為數值變數建立字串版本的變數
  • 利用過程步中的選項

（1）ROUND函式
在不需要特別精確的情況下，可以使用ROUND函式僅比較前幾位小數。

data s3;
    merge s1 s2;
    by index;
    a_r = round(a, 0.1);
    b_r = round(b, 0.1);
    if a_r = b_r then A_eq_B = "Y";
run;
proc print data = s3;
    var index a b A_eq_B a_r b_r;
    format a hex16. b hex16. a_r hex16. b_r hex16.;
run;

（2）建立數值變數的字串版本

data s3;
    merge s1 s2;
    by index;
    a_fmt = put(a, 3.1);
    b_fmt = put(b, 3.1);
    if a_fmt = b_fmt then A_eq_B = "Y";
run;
proc print data = s3;
    var index a b A_eq_B a_fmt b_fmt;
    format a hex16. b hex16. a_fmt $hex16. b_fmt $hex16.;
run;

（3）使用過程步中的選項

proc compare base = s1 compare = s2(rename = (b = a)) criterion = 0.1 method = absolute;
run;

參考文獻：https://www.jianguoyun.com/p/DfN3qOQQ6YbXCRimmbME

相關文章