統計学対数変換後、元のスケールに戻して信頼区間を算出する方法

Question

IPAが公開しているソフトウエア開発データ白書で、
信頼区間　５０％、９５％の具体的な計算方法がわかりません。
計算で信頼区間を算出し、散布図に描画したのですが、公開されているデータ、図とぜんぜん合いません。


対象データ：
https://www.ipa.go.jp/ikc/publish/whitepaper_dl.html

対象ファイル：
graph_data_06.xlsx
「6-2-1」シート

対応する図：
ソフトウェア開発データ白書 2018-2019_本編.pdf
P136 図表6-2-2


■やったこと：
信頼区間の算出については、以下のサイトを参考にしました。
https://cyclo-commuter.hatenablog.jp/entry/2018/02/02/090255

エクセルで、データ白書のデータについて対数変換してから、
信頼区間をプロットしてみたのですが、IPAの図表の信頼区間と全く合致しません。。
（データと近似曲線（累乗近似）は、エクセルでIPAどおりに描画できています）

統計については、エクセルやR言語、EZRで最近勉強はじめまた初心者です。
エクセル、R言語、どちらの方法でも結構ですので、対数スケールに変換してから通常スケールに戻す場合の信頼区間の計算方法についてご教示ください。

■補足１
X軸:実績工数（人時）
Y軸:実績月数

データ白書の説明（p28-29）によれば、プロジェクトのデータは正規分布していない事が多いが、
対数に変換するとほぼ正規分布と見なせる。このため、対数スケールに変換してから通常スケールに
戻して信頼区間付き散布図を作成したとあります。

■補足２
データ白書のP29に「※ 信頼区間の算出は、
田中豊、脇本和昌、｢多変量統計解析法｣、現代数学社、1983 
涌井良幸、涌井貞美、｢図解でわかる回帰分析、日本実業出版社、2002 を参考にした。」
とあり、まず「図解でわかる回帰分析」を読んでいるところですが、道は長いです。
非線形回帰の話が出てくるのですが、単回帰で計算できるはずなので、深いところまでは
追っていません。統計は奥が深いですね。

kamiyasiro · Accepted Answer

#1です。

ちなみに、yを一旦対数変換し、等分散仮定で、ｙ～log(x)のモデルで近似式を作り、全体をexpで戻してプロットした図を添付します。

ご質問者様も、こんなグラフになったのではないでしょうか？

kamiyasiro · Answer

#1です。

信頼限界幅の係数ですが、参考にされたＨＰの式は係数がルートの中に入っているので、F値で良いです。

kamiyasiro · Answer

#1です。

rm(list = ls())

x <- read.csv("6-2-1.csv")
colnames(x) <- c("実績工数（プロジェクト全体）［人時］",
                 "実績月数（プロジェクト全体）［月］")

par(mar=c(10.1,4.1,8.1,2.1))
plot(x, pch = 4, col = 4, cex = 0.6, cex.lab = 0.8, las = 1,
     xaxs = "i", yaxs = "i", xlim = c(0, 300000), ylim = c(0, 40),
     main = "図表 6-2-1 ● プロジェクト全体の工数と工期（新規開発）（信頼区間 50%、95% 付き）",
     cex.main = 0.8)
axis(side = 2, tck = 1.0, labels = FALSE)

n <- nrow(x)
mtext(paste("N =",n), line = 0, at = 280000)
phi <- n - 2    # 原点を通しているから１かもしれない

# -------------------------------------------------------
# GLM

colnames(x) <- c("x", "y")
result <- glm(y ~ log(x), data = x, family = "Gamma")
summary(result)

# 近似線を求める計算点を準備
xpred <- data.frame(x = seq(0, 300000, by=1000))

# 近似値の計算と描画

d <- 1/summary(result)$dispersion
d

library("MASS")
d <- gamma.shape(result)$alpha
d      # こちらの方が良いらしい

preds <- predict(result, xpred, type = "response", dispersion = d, se.fit=TRUE)
lines(cbind(xpred, preds$fit))

# 95%信頼区間
critp <- qt(0.975, phi)
lines(cbind(x = xpred, y = preds$fit + critp * preds$se.fit), col = 2)
lines(cbind(x = xpred, y = preds$fit - critp * preds$se.fit), col = 2)

# 50%信頼区間
critp <- qt(0.75, phi)
lines(cbind(x = xpred, y = preds$fit + critp * preds$se.fit), col = 6)
lines(cbind(x = xpred, y = preds$fit - critp * preds$se.fit), col = 6)

# 凡例を描く（省略）

par(mar=c(5.1,4.1,4.1,2.1))                   #デフォルト

kamiyasiro · Answer

#1です。

対数を取ってから平均を取るために足し算すると、対数の足し算は掛け算だから幾何平均を取っていることになるのです。

回帰に関してはガウスマルコフの定理っていうのがあり、それを満たさなければなりません。そのひとつにばらつきの期待値E(ε)は０、つまり上下の平均を通るというのがあるのですが、幾何平均を通るとガウスマルコフの定理に反し、不偏推定量にならないのです。

あの図6-2-1の信頼区間は絶対におかしいです。あれが再現できても意味ないですよ。いいですか。信頼区間っていうのは、いわば回帰線の存在範囲ですよ。にもかかわらず、あの図の95％信頼区間は、まるで予測区間（データの存在範囲）のような位置に来ています。

手順に則って描くと添付図のようになります。誤差モデルはGamma（変動係数一定）にしています。

続いて、Rスクリプトを投稿します。

kamiyasiro · Answer

#1です。

普通はyがリニアスケールのときは、回帰線は平均を通っていきますよね。

ｙを対数変換したときに平均を通るってことは、expで戻したときは幾何平均を通るってことですよね。そんなの回帰線の定義から外れますよね。

だから、あのグラフ、回帰線をあえて入れてないのですかね。

kamiyasiro · Answer

#1です。

『データ白書の説明（p28-29）によれば、プロジェクトのデータは正規分布していない事が多いが、対数に変換するとほぼ正規分布と見なせる。このため、対数スケールに変換してから通常スケールに戻して信頼区間付き散布図を作成したとあります。』

これって、ｙの誤差の話？

だったら、上下非対称というのも頷けます。

ｙを対数変換した空間で回帰し信頼区間を求めたうえで、expで全体を戻すんですかね。
そのときも、ｙ～A＋x^B の式で良いのかしらん？

誤差モデルとしては、なんだか変な話ですよね。

面白い問題ですね。

kamiyasiro · Answer

#1です。

いやー、あの図の信頼区間はおかしいですよ。なぜなら、ｙはリニアスケールなのに、上下で非対称！

poissonでやってもGAMMAでやっても合わないので、2時間くらい悩んで無駄な時間を費やしてしまいました。

なぜ、最初に気づかなかったのだろう。ちゃんとした解析者がやったと思って結果を信じたから？

ご質問者様が合わないのも無理ないですよ。

あとで私が作成したRスクリプトを提供します。

kamiyasiro · Answer

企業で統計を推進する立場の者です。

グラフ見ました。

一般化線形モデルでやっているとみて、今、Rスクリプト書き始めたんですが、プロットが違います。対応する図は6-2-2ではなく6-2-1ですね。

ご質問様のやられている数値変換して線形回帰するのとは全く異なる方法です。誤差の仮定が違います。

一致する図ができたら、Rスクリプトを提供します。

統計学 対数変換後、元のスケールに戻して信頼区間を算出する方法

#1です。

#1です。

#1です。

#1です。

#1です。

#1です。

#1です。

企業で統計を推進する立場の者です。

似たような質問が見つかりました

関連するカテゴリからQ&Aを探す

デイリーランキングこのカテゴリの人気デイリーQ&Aランキング

マンスリーランキングこのカテゴリの人気マンスリーQ&Aランキング

統計学対数変換後、元のスケールに戻して信頼区間を算出する方法