メイン mreg-qt1.html Last modified: Sep 01, 2009
<html>
<head>
<meta http-equiv="Content-Type" content="text/html;CHARSET=EUC-JP">
<link rel="shortcut icon" href="../favicon.ico">
<title>JavaScript</title>
<script src="inv.js">document.write("inv.js ファイルが見つかりません??<br>")</script>
<script src="fxp.js">document.write("fxp.js ファイルが見つかりません??<br>")</script>
<script src="io.js">document.write("io.js ファイルが見つかりません??<br>")</script>
<script language="JavaScript">
<!--
var epsinv = 1e-8
function zfi(w, i)
{
var ww = (""+i).length
ww = ww > w ? 0 : w-ww
return "000000000000000".substring(0, ww)+i
}
function vname(i)
{
return "V"+zfi(2, i+1)
}
function sort(x, n, ip)
{
var min, minp, i, j, temp
for (i = 0; i < n-1; i++) {
min = x[i][ip]
minp = i
for (j = i+1; j < n; j++) {
if (x[j][ip] < min) {
min = x[j][ip]
minp = j
}
}
if (i != minp) {
temp = x[i][ip]
x[i][ip] = min
x[minp][ip] = temp
}
}
for (i = 0; i < n; i++) {
if (x[i][ip] != i+1) {
return 0
}
}
return 1
}
function itemno(data, nc, item, cat, junjo, value)
{
var i, j, k, nobe, flag0 = 0, flag, work, catj
for (j = 0; j < item; j++) {
catj = flag = 0
for (i = 0; i < nc; i++) {
work = data[i][j]
if (Math.floor(work) != work) {
if (flag == 0) {
printf("アイテム変数 %s が，整数値でありません。\n", vname(j))
flag = flag0 = 1
}
}
value[catj][j] = work
for (k = 0; k <= catj; k++) {
if (value[k][j] == work) {
break
}
}
if (k == catj) {
catj++
if (catj > 100) {
printf("アイテム変数 %s の取る値の種類が多すぎます。\n", vname(j))
flag0 = 1
break
}
}
}
cat[j] = catj
if (catj == 1) {
printf("アイテム変数 %s は，一種類の値しか取りません。この変数は分析に使用しないでください。\n", vname(j))
flag0 = 1
}
else if (sort(value, catj, j) == 0) {
for (i = 0; i < nc; i++) {
work = data[i][j]
for (k = 0; k < catj; k++) {
if (work != k+1 && work == value[k][j]) {
data[i][j] = k+1
break
}
}
}
}
}
if (flag0) {
return -99
}
else {
nobe = 0
for (i = 0; i < item; i++) {
junjo[i] = nobe
nobe += cat[i]
}
return nobe
}
}
function outstat(mean, variance, sd, r, nv, names)
{
var j, k
printf("%10s%15s%15s%15s\n", "変数", "平均値", "不偏分散", "標準偏差")
for (j = 0; j < nv; j++) {
printf("%10s%15.7g%15.7g%15.7g\n", names[j], mean[j], variance[j], sd[j])
}
printf("\n相関係数行列\n\n ")
for (j = 0; j < nv; j++) {
printf("%8s", names[j])
}
printf("\n")
for (j = 0; j < nv; j++) {
printf("%10s", names[j])
for (k = 0; k < nv; k++) {
printf("%8.3f", r[j][k])
}
printf("\n")
}
}
function getdummy(i, data, item, nobe, junjo, dd)
{
var j, k
for (j = 0; j < nobe; j++) {
dd[j] = 0
}
for (j = 0; j < item; j++) {
k = junjo[j]+data[i][j]-1
dd[k] = 1
}
dd[nobe] = data[i][item]
}
function mreg(flag, nc, nv, mean, variance, sd, r, ss, names, freq, data, item, nobe, junjo, order, order2, use, dd, value, method)
{
var b, b1, stde, nv1, nc1, b0, se, i, j, temp2, sr, rhoy, ve, retv, adj
var syy, seb0, df2, seb, p, ip, pred, stres, j1, sum, sum2
nv1 = nv-1
nc1 = nc-1
if (flag) {
if (nv == 1) {
printf("基準を満たす独立変数がありません\n")
return
}
printf("★ 最もよい重回帰モデル ★\n\n")
printf("有効データ組数 = %i\n\n", nc)
printf("独立変数の個数 = %i\n\n", nv1)
outstat(mean, variance, sd, r, nv, names)
}
b = new Array(nv)
b1 = new Array(nv)
stde = new Array(nv)
if (inv(r, nv1, epsinv) == 0) {
b0 = mean[nv1]
se = syy = ss[nv1]
for (i = 0; i < nv1; i++) {
temp2 = 0
for (j = 0; j < nv1; j++) {
temp2 += r[i][j]*r[j][nv1]
}
b1[i] = temp2
b[i] = temp2*sd[nv1]/sd[i]
b0 -= b[i]*mean[i]
se -= b[i]*ss[i]
}
sr = syy-se
rhoy = se/syy
seb0 = 0
df2 = nc-nv1-1
for (i = 0; i < nv1; i++) {
seb = r[i][i]*rhoy/df2*(Math.pow(sd[nv1]/sd[i],2))
seb0 += mean[i]*mean[i]*seb
stde[i] = Math.sqrt(seb)
}
seb = 0
for (i = 0; i < nv1; i++) {
temp2 = 2*mean[i]/sd[i]
for (j = i+1; j < nv1; j++) {
seb += r[j][i]*temp2*mean[j]/sd[j]
}
}
seb0 = Math.sqrt(seb0+(seb+nc1/nc)*sd[nv1]*sd[nv1]*rhoy/df2)
if (flag) {
printf("\n%10s%15s%15s%12s%12s%18s\n", "変数", "偏回帰係数", "標準誤差", "t値", "P値", "標準化偏回帰係数")
for (i = 0; i < nv1; i++) {
if (stde[i] == 0) {
printf("偏回帰係数の標準誤差がゼロになります（従属変数が完全に予測できます）。このようなことは普通はありません。特殊なデータですね。意図して入力したデータですか？\n")
return
}
p = fxp(Math.pow(b[i]/stde[i], 2), 1, df2)
printf("%10s%15.7g%15.7g%12.5f%12.5f%15.7g\n", names[i], b[i], stde[i], Math.abs(b[i])/stde[i], p, b1[i])
}
printf("%10s%15.7g%15.7g%12.5f%12.5f\n", "定数項", b0, seb0, b0/seb0, fxp(Math.pow(b0/seb0, 2), 1, df2))
printf(" %s = %i\n\n", "t値の自由度", df2)
ve = se/df2
temp2 = sr*df2/nv1/se
p = fxp(temp2, nv1, df2)
printf("分散分析表\n\n")
printf("%4s%15s%10s%15s%12s%12s\n", "要因", "平方和", "自由度", "平均平方", "F値", "P値")
printf("%4s%15.7g%10i%15.7g%12.5f%12.5f\n", "回帰", sr, nv1, sr/nv1, temp2, p)
printf("%4s%15.7g%10i%15.7g\n", "残差", se, df2, ve)
printf("%4s%15.7g%10i%15.7g\n\n", "全体", syy, nc1, syy/nc1)
}
adj = 1-se/df2/syy*nc1
retv = (method == 1) ? 1-rhoy : adj
if (flag) {
printf("重相関係数 = %.5f 決定係数 = %.5f 修正決定係数 = %.5f\n\n", Math.sqrt(1-rhoy), 1-rhoy, adj)
}
if (flag) {
printf("数量化 I 類の出力形式による結果\n\n")
printf("%20s %19s\n", "アイテム・カテゴリー", "カテゴリースコア")
sum = new Array(item)
for (i = 0; i < item; i++) {
sum[i] = 0
}
sum2 = 0
for (i = 0; i < nv1; i++) {
sum[order2[i]] += freq[i]*b[i]
}
ip = 0
for (i = 0; i < nobe; i++) {
if (order[i] == -1 && use[order[i+1]] == 1) {
j = order[i+1]
printf("%14s %21.7g\n", vname(order[i+1])+"/"+value[0][j], -sum[j]/nc)
sum2 += sum[j]/nc
}
else if (order[i] != -1 && use[order[i]] == 1) {
j = order[i]
printf("%14s %21.7g\n", names[ip], b[ip]-sum[j]/nc)
ip++
}
}
printf("%14s %21.7g\n\n", "定数項", b0+sum2)
}
if (flag) {
for (j = 0; j < nv-1; j++) {
sd[j] *= Math.sqrt((nc-1)/nc)
}
printf("予測結果\n\n%5s %12s %12s %12s %12s\n", "番号", "観察値", "予測値", "残差", "標準化残差")
for (i = 0; i < nc; i++) {
getdummy(i, data, item, nobe, junjo, dd)
ip = 0
for (j = 0; j <= nobe; j++) {
if (j == nobe || (order[j] != -1 && use[order[j]] == 1)) {
dd[ip++] = dd[j]
}
}
pred = b0
for (j = 0; j < nv-1; j++) {
pred += b[j]*dd[j]
b1[j] = (dd[j]-mean[j])/sd[j]
}
stres = 0
for (j = 0; j < nv-1; j++) {
for (j1 = 0; j1 < nv-1; j1++) {
stres += b1[j]*b1[j1]*r[j1][j]
}
}
stres = (dd[j]-pred)/Math.sqrt(ve*(1-(stres+1)/nc))
printf("%5i %12.6f %12.6f %12.6f %12.6f\n", i+1, dd[j], pred, dd[j]-pred, stres)
}
}
}
return retv
}
function submatrix(nobe, use2, order, order2, mean, mean2, variance, variance2, sd, sd2, ss, ss2, r, r2, names, names2, freq, freq2)
{
var ip, jp, i, j
ip = 0
for (i = 0; i <= nobe; i++) {
if (i == nobe || use2[order[i]] == 1) {
order2[ip] = order[i]
freq2[ip] = freq[i]
mean2[ip] = mean[i]
variance2[ip] = variance[i]
sd2[ip] = sd[i]
ss2[ip] = ss[i]
names2[ip] = names[i]
jp = 0
for (j = 0; j <= nobe; j++) {
if (j == nobe || use2[order[j]] == 1) {
r2[ip][jp] = r[i][j]
jp++
}
}
ip++
}
}
return ip
}
function calc(data_string, method, alpha0)
{
var data, nc, nv, i, j, k, mean, variance, sd, ss, r, mean2, variance2, sd2, ss2, r2, dd, cat, use, use2
var item, junjo, value, nobe, order, temp, loop, loop0, step, freq, freq2, names2, order2, used, nused
var r_sq, step, r_sqdif, fdif, pval
var ip, stat, maxp, max, oldmax, names
if (isNaN(alpha = parseFloat(alpha0)) || alpha <= 0 || alpha >= 1) {
alpha = 0.05
}
if ((data = getdata(data_string, 0)) == false) return
nc = data.length
nv = data[0].length
if (nc < 2) {
printf("ケース数が１以下です\n")
return
}
else if (nc <= nv) {
printf("データの組数は，独立変数の個数より2以上大きくないといけません\n")
return
}
item = nv-1
cat = makeVector(item)
junjo = makeVector(item)
value = makeMatrix(100, item)
item = nv-1
nobe = itemno(data, nc, item, cat, junjo, value)
if (nobe <= 0) {
printf("データをチェックしてください\n")
return
}
else if (nobe >= 500) {
printf("カテゴリーが多すぎます\n")
return
}
freq = new Array(nobe+1)
mean = new Array(nobe+1)
variance = new Array(nobe+1)
r = makeMatrix(nobe+1, nobe+1)
sd = new Array(nobe+1)
ss = new Array(nobe+1)
names = new Array(nobe+1)
order = new Array(nobe+1)
freq2 = new Array(nobe+1)
mean2 = new Array(nobe+1)
variance2 = new Array(nobe+1)
r2 = makeMatrix(nobe+1, nobe+1)
sd2 = new Array(nobe+1)
ss2 = new Array(nobe+1)
names2 = new Array(nobe+1)
order2 = new Array(nobe+1)
dd = new Array(nobe+1)
use = new Array(item)
use2 = new Array(item)
used = new Array(item)
// Netscape communicator のバグ？配列ベクトルは初期化されていない
for (i = 0; i <= nobe; i++) {
freq[i] = 0
mean[i] = 0
for (j = 0; j <= i; j++) {
r[i][j] = 0
}
}
// ダミー変数と外的基準の平均値・分散・標準偏差と相関係数行列
for (i = 0; i < nc; i++) {
getdummy(i, data, item, nobe, junjo, dd)
for (j = 0; j <= nobe; j++) {
if (dd[j] == 1) {
freq[j]++
}
mean[j] += dd[j]
}
}
for (j = 0; j <= nobe; j++) {
mean[j] /= nc
}
for (i = 0; i < nc; i++) {
getdummy(i, data, item, nobe, junjo, dd)
for (j = 0; j <= nobe; j++) {
temp = dd[j]-mean[j]
for (k = 0; k <= j; k++) {
r[j][k] += temp*(dd[k]-mean[k])
}
}
}
for (j = 0; j <= nobe; j++) {
ss[j] = r[nobe][j]
variance[j] = r[j][j]/(nc-1)
sd[j] = Math.sqrt(variance[j])
for (k = 0; k <= j; k++) {
r[j][k] /= (nc-1)*sd[j]*sd[k]
r[k][j] = r[j][k]
}
}
// ダミー変数の所属元
for (i = 0; i < item; i++) {
order[junjo[i]] = -1
names[junjo[i]] = vname(i)+"/"+value[0][i]
for (j = 1; j < cat[i]; j++) {
order[junjo[i]+j] = i
names[junjo[i]+j] = vname(i)+"/"+value[j][i]
}
}
order[nobe] = item
names[nobe] = vname(item)
printf("★★ カテゴリー変数の選択（変数増加法）を行う重回帰分析（数量化 I 類） ★★\n\n")
for (i = 0; i < item; i++) {
use[i] = 0
}
max = 0
nused = 0
r_sq = 0
step = 0
for (loop0 = 0; loop0 < item; loop0++) {
oldmax = max
for (loop = 0; loop < item; loop++) {
for (i = 0; i < item; i++) {
use2[i] = use[i]
}
if (use2[loop] == 0) {
use2[loop] = 1
ip = submatrix(nobe, use2, order, order2, mean, mean2, variance, variance2, sd, sd2, ss, ss2, r, r2, names, names2, freq, freq2)
stat = mreg(0, nc, ip, mean2, variance2, sd2, r2, ss2, names2, freq2, data, item, nobe, junjo, order, order2, use2, dd, value, method)
if (stat > max) {
max = stat
maxp = loop
}
}
}
if (method == 1) {
r_sqdif = max-r_sq
fdif = r_sqdif*(nc-ip)/(cat[maxp]-1)/(1-max)
pval = fxp(fdif, cat[maxp]-1, nc-ip)
if (pval > alpha) {
break
}
r_sq = max
}
else if (max == oldmax) {
break
}
used[nused++] = maxp
use[maxp] = 1
printf("** 第 %i ステップ **\n", ++step)
printf("追加された変数 = %s\n", vname(maxp))
printf("使用された変数 : ")
for (i = 0; i < nused-1; i++) {
printf("%s, ", vname(used[i]))
}
printf("%s\n", vname(used[i]))
if (method == 1) {
printf("決定係数 = %g\n 決定係数の増加分 = %g\n 決定係数の増加分の F 値 = %g 自由度 = (%i, %i)\n 決定係数の増加分の P 値 = %g （α = %g）\n\n", max, r_sqdif, fdif, cat[maxp]-1, nc-ip, pval, alpha)
}
else {
printf("修正決定係数 = %g\n\n", max)
}
}
ip = submatrix(nobe, use, order, order2, mean, mean2, variance, variance2, sd, sd2, ss, ss2, r, r2, names, names2, freq, freq2)
mreg(1, nc, ip, mean2, variance2, sd2, r2, ss2, names2, freq2, data, item, nobe, junjo, order, order2, use, dd, value, method)
}
//-->
</script>
</head>
<body bgcolor="#ffffff">
<font size="+2"><b>カテゴリー変数の選択（変数増加法）を行う重回帰分析（数量化 I 類）</b></font> <a href="src/mreg-qt1.html"><img src="png/src.png" width=35 height=11 alt="src" align=top></a> Last modified: Jun 01, 2006<hr noshade><p>
<font color="#ff0000" size="+2">以下のプログラムのサポートは終了しました。自己責任でお使い下さい。</font>
<form name=Result>
<script language="JavaScript">
<!--
//-->JavaScript がサポートされていないブラウザですか？
</script>
<table><tr>
<td><input type="button" name="calcurate" value="計算開始" onClick="calc(this.form.data.value, this.form.method.value, this.form.alpha.value)"> </td>
<td><input type="button" name="clear" value="入力欄クリア" onClick="this.form.data.value=''"> </td>
<td><input type="button" name="clear" value="出力欄クリア" onClick="this.form.result.value=''"> </td>
<td nowrap><img src="../gra/button3.png" width=9 height=9 alt="・"> <a href="exa/mreg-qt1.html">使用法</a></td>
</tr></table>
<p>
解の探索方法<input name="method" value="1" size=5>（1: 決定係数の増分の有意性， 2: 修正決定係数が増加する限り）<br>
有意水準（1の場合）<input name="alpha" value="0.05" size=5><br>
入力欄<br><textarea name="data" ROWS=20 COLS=80></textarea><p>
出力欄<br><textarea name="result" ROWS=30 COLS=80></textarea>
</form>
<p><hr noshade>
<img src="../gra/button3.png" width=9 height=9 alt="・"> <a href="../lecture/Regression/mreg/index.html">手法の解説(1)</a> <img src="../gra/button3.png" width=9 height=9 alt="・"> <a href="../lecture/Qt/qt1.html">手法の解説(2)</a><br>
<img src="../gra/button3.png" width=9 height=9 alt="・"> <A HREF="javascript:history.go(-1)">直前のページへ戻る</A> <img src="../gra/button3.png" width=9 height=9 alt="・"> <a href="../mail.html">E-mail to Shigenobu AOKI</a>
<p><center><IMG SRC="../gra/ume5.png" width=121 height=37 ALT="Made with Macintosh"></center>
</body>
</html>
サブ inv.js Last modified: Mar 25, 2004
サブ fxp.js Last modified: Mar 25, 2004
サブ io.js Last modified: Mar 25, 2004
直前のページへ戻る E-mail to Shigenobu AOKI
