// 以下是计算t分布及线性回归等程序 function normsinv(pp){ var x,y; var q = 1 - pp; bb = [1.570796288,0.03706987906,-0.0008364353589,-0.0002250947176, 0.000006841218299,0.000005824238515,-0.00000104527497,8.360937017e-8,-3.231081277e-09, 3.657763036e-11,6.936233982e-13]; if (q == 0.5) return 0; var p; if (q > 0.5) p = 1 - q; else p = q; var b; y =-Math.log(4 * p * (1 - p)); b = y * bb[10]; var i; for (i = 9; i >= 0; i--) { b = y * (bb[i] + b); } x = Math.sqrt(b); if (q > 0.5) x = -x; return x; } function gamma( x) { var i; var y, t, s, u; var a = [ 0.0000677106,-0.0003442342, 0.0015397681,-0.0024467480,0.0109736958, -0.0002109075,0.0742379071,0.0815782188, 0.4118402518,0.4227843370,1.0]; if (x <= 0.0) throw "x不能小于0"; y = x; if (y <= 1.0) { t = 1.0 / (y * (y + 1.0)); y += 2.0; } else if (y <= 2.0) { t = 1.0 / y; y += 1.0; } else if (y <= 3.0) t = 1.0; else { t = 1.0; while (y > 3.0) { y -= 1.0; t = t * y; } } s = a[0]; u = y - 2.0; for (i = 1; i <= 10; i++) s = s * u + a[i]; s *= t; return s; } function betacf(a, b, x) { var MAXIT = 1000; var EPS = 1e-10; var FPMIN = 1e-8; var m, m2; var aa, c, d, del, h, qab, qam, qap; qab = a + b; qap = a + 1; qam = a - 1; c = 1; d = 1 - qab * x / qap; if (Math.abs(d) < FPMIN) d = FPMIN; d = 1 / d; h = d; for (m = 1; m <= MAXIT; m++) { m2 = 2 * m; aa = m * (b - m) * x / ((qam + m2) * (a + m2)); d = 1 + aa * d; if (Math.abs(d) < FPMIN) d = FPMIN; c = 1 + aa / c; if (Math.abs(c) < FPMIN) c = FPMIN; d = 1 / d; h *= d * c; aa = -(a + m) * (qab + m) * x / ((a + m2) * (qap + m2)); d = 1 + aa * d; if (Math.abs(d) < FPMIN) d = FPMIN; c = 1 + aa / c; if (Math.abs(c) < FPMIN) c = FPMIN; d = 1 / d; del = d * c; h *= del; if (Math.abs(del - 1) <= EPS) break; } if (m > MAXIT) throw "a or b too big, or MAXIT too small in betacf"; return h; } function gammln(xx) { var j; var x, y, tmp, ser; var cof = [76.18009172947146,-86.50532032941677,24.01409824083091,-1.231739572450155, 0.1208650973866179e-2,-0.5395239384953e-5]; y = x = xx; tmp = x + 5.5; tmp -= (x + 0.5) * Math.log(tmp); ser = 1.000000000190015; for (j = 0; j < 6; j++) ser += cof[j] / ++y; return -tmp + Math.log(2.5066282746310005 * ser / x); } function betai(a, b, x) { var bt; if (x < 0 || x > 1) throw "Bad x in routine batai"; if (x == 0 || x == 1) bt = 0; else bt = Math.exp(gammln(a + b) - gammln(a) - gammln(b) + a * Math.log(x) + b * Math.log(1 - x)); if (x < (a + 1) / (a + b + 2)) return bt * betacf(a, b, x) / a; else return 1 - bt * betacf(b, a, 1 - x) / b; } function tdist2(t, free) { if (free < 1) throw "free must biger than 0"; var f; f = free; var p; p = 1 - betai(f / 2.0, 0.5, f / (f + t * t)); return p; } function tdist(t, free) { if (t == 0) { return 0.5; } if (t > 0) { return 0.5 + tdist2(t, free) / 2; } return (1 - tdist2(-t, free)) / 2; } var outpp; function t_dist(n, t) { var sign, ibi, n2, i; var tt, x, p, u, ga1, ga2,dd; if (t == 0) { ga1 = gamma(n / 2.0); ga2 = gamma(0.5 + n / 2.0); outpp = 0.5; return ga2 / (Math.sqrt(n * Math.PI) * ga1); } if (t < 0) { sign = -1; } else sign = 1; tt = t * t; x = tt / (n + tt); if (n % 2 == 0) { p = Math.sqrt(x); u = p * (1 - x) / 2.0; ibi = 2; } else { u = Math.sqrt(x * (1 - x)) / Math.PI; p = 1 - 2 * Math.atan(Math.sqrt((1 - x) / x)) / Math.PI; ibi = 1; } if (ibi == n) { dd = u / Math.abs(t); outpp = 0.5 + sign * p / 2.0; return dd; } else { n2 = n - 2; } for (i = ibi; i <= n2; i += 2) { p += 2.0 * u / i; u = u * (1 + i) / i * (1 - x); } dd = u / Math.abs(t); outpp = 0.5 + sign * p / 2.0; return dd; } function tdistinv(pvar, n) { var q; q = 1 - pvar; var t=0.0; var pis, dfr2, c, q2, p, yq, e, ga1, ga2, ga3, t0, pp, d; var k; pis = Math.sqrt(Math.PI); dfr2 = n / 2.0; if (n == 1) return Math.tan(Math.PI * (0.5 - q)); if (n == 2) { if (q > 0.5) c = -1; else c = 1; q2 = (1 - 2 * q); q2 *= q2; return Math.sqrt(2.0 * q2 / (1.0 - q2)) * c; } p = 1 - q; yq = normsinv(p); e = 1.0 - 1.0 / (4.0 * n); e = e * e - yq * yq / (2.0 * n); if (e > 0.5) t0 = yq / Math.sqrt(e); else { ga1 = gammln(dfr2); ga2 = gammln(dfr2 + 0.5); ga3 = Math.exp((ga1 - ga2) / n); t0 = Math.pow(Math.sqrt(n) / (pis * q * n), (1.0 / n) / ga3); } for (k = 1; k <= 30; k++) { d = t_dist(n, t0); pp=outpp; if (d == 0) return t0; t = t0 - (pp - p) / d; if (Math.abs(t0 - t) < 0.000001 * Math.abs(t)) { return t; } else t0 = t; } return t; } function relhold(freedom, alpha) { var t; t = tdistinv(1 - alpha / 2, freedom); return t / Math.sqrt(freedom + t * t); } // 线性回归的核心程序 function linehuigui(data, alpha) { var a,b,sig2,sigx2,r,ra,t,ta,xmean; var n; n = data.length / 2; var i; var sumx = 0, sumy = 0, sumxx = 0, sumyy = 0, sumxy = 0; var x, y; for (i = 0; i < n; i++) { x = data[2 * i]; y = data[2 * i + 1]; sumx += x; sumy += y; sumxx += x * x; sumyy += y * y; sumxy += x * y; } xmean = sumx / n; var sxx, syy, sxy; sxx = sumxx - sumx * sumx / n; syy = sumyy - sumy * sumy / n; sxy = sumxy - sumx * sumy / n; b = sxy / sxx; a = (sumy / n) - (sumx / n) * b; sig2 = (syy - b * sxy) / (n - 2.0); sigx2 = sxx; t = b * Math.sqrt(sxx / sig2); ta = tdistinv(1 - alpha / 2, n - 2); r = sxy / Math.sqrt(sxx * syy); ra = relhold(n - 2, alpha); var out=new Array(); out[0]=a; out[1]=b; out[2]=sig2; out[3]=sigx2; out[4]=r; out[5]=ra; out[6]=t; out[7]=ta; out[8]=xmean; return out; } // 以下是绘制曲线程序 function updown(a){// 将给定的数分离成两个数相乘，一个是数的部分，另一个是指数部分，利用函数toExponential来做 var r=new Array(2); var s,s1,s2; s=a.toExponential(4); var i,j,c,b; i=s.indexOf("e"); j=s.length; s1=s.substring(0,i); s2=s.substring(i+1,j); c=parseFloat(s1); b=parseFloat(s2); r[0]=c; r[1]=b; return r; } function getscals(low, high)// 根据最低值low和最高值high计算刻度的数量和刻度，返回一个实数矩阵，维数就是刻度，内容就是各个刻度的值 { var l=high-low; if(low>0)l=high; if(high<0){ l=-low; high=0; } var llow=low; var lhigh=high; var bs=new Array(); var bbb=new Array(); var bbb=updown(l/4); var a=bbb[0]; var index=bbb[1]; llow/=Math.pow(10,index); lhigh/=Math.pow(10,index); var ascale; if(a<=2) ascale=2; else if(a<=5){ ascale=5; } else ascale=10; var abegin=0; var i,d,dd,s; if(llow<0){ for(i=1;i<1000;i++){ abegin=-i*ascale; if(abegin<=llow) break; } } else abegin=0; for(i=0;i<50;i++){ d=abegin+i*ascale; dd=d*Math.pow(10,index); s=dd.toPrecision(8); bs[i]=parseFloat(s); //bs[i]=d*Math.pow(10,index); if(d>=lhigh || Math.abs(d-lhigh)<1e-5)break; } return bs; } var ctx; var width=600,hight=400;//图像的宽度和高度 var xbegin,ybegin,xlength,ylength;//曲线的绘制范围，x和y的起始点和宽高度 var axmin,aymin,axwidth,aywidth; function drawaxis(xl,xh,yl,yh){//绘制xy坐标轴，xl,xh,yl,yh是x和y取值的最小值和最大值 var yofx,xofy;//,ybegin,xbegin;//x轴的y坐标和y轴的x坐标及曲线的x和y的起始点 var xscalelength,yscalelength;//两个刻度之间的距离 //var xlength,ylength;//曲线的x和y的幅度 var xscals,yscals; // 是数组，存放x和y的刻度 var iyofx,ixofy;//y刻度的第几个数是0，x刻度的第几个数是0 xscals=getscals(xl,xh); axmin=xscals[0]; axwidth=xscals[xscals.length-1]-xscals[0]; yscals=getscals(yl,yh); aymin=yscals[0]; aywidth=yscals[yscals.length-1]-yscals[0]; yofx=ybegin=hight-30; ylength=hight-70-30; yscalelength=ylength/(yscals.length-1); var i; for(i=0;ij;i--){ // 删去末尾的字符0 if(s.substr(i,1)!="0"){ return s.substr(0,i+1); } } return s.substr(0,j); } function strtodouble(s) { var index=s.search(/\//); if(index<0){ return parseFloat(s); } var ss="abc"; ss=s.substring(index+1,s.length); var a,b; var sss; sss=s.substring(0,index); a=parseFloat(ss); b=parseFloat(sss); //return (0.0+sss)/(0.0+ss); return b/a; } function stringstodouble(s) { var reg_pattern=/[0123456789.\-\/]+/g; var arr=s.match(reg_pattern); var a=new Array(); var i; if(arr!=null) for(i=0;ixmax)xmax=datax[i]; if(datay[i]ymax)ymax=datay[i]; } var canvas=document.getElementById('my_canvas'); if(canvas.getContext) { ctx=canvas.getContext('2d'); ctx.globalCompositeOperation="source-over"; ctx.clearRect(0,0,600,400); ctx.lineWidth=1; // 在下面开始绘制曲线 drawaxis(xmin,xmax,ymin,ymax); ctx.closePath(); ctx.beginPath(); var ix,iy; for(i=0;i"; spanresult.innerHTML+="从t检验的角度看：t检验统计量为t=" +dtos(t) + ",而t_{" +dtos(alpha/2) + "}(" + (n - 2).toString() + ")=" +dtos(ta)+"
"; if (Math.abs(t) > ta) { spanresult.innerHTML += "因为|t|>t_α/2, 所以拒绝b=0的原假设，回归效果显著。"; divpre.style.display="block"; } else { spanresult.innerHTML += "因为|t|≤t_α/2, 所以接受b=0的原假设，回归效果不好。"; } spanresult.innerHTML+= "
从相关系数检验法的角度看：r=" + dtos(r) + ",r_{" +dtos(alpha) + "}=" + dtos(ra)+"
"; if (Math.abs(r) > ra) { spanresult.innerHTML += "因为|r|>r_α, 所以拒绝b=0的原假设，回归效果显著。"; } else { spanresult.innerHTML += "因为|r|≤r_α, 所以接受b=0的原假设，回归效果不好。"; } ctx.lineWidth=1; ctx.beginPath();//注意，从此开始的四句毫无意义，但是如果不加，就会有重叠 ctx.moveTo(1,1); ctx.lineTo(2,1); ctx.closePath(); ctx.stroke(); } } function Button2_onclick() { var yinc, x0, y0; if(Text1.value=="" || Text1.value==null){ spanr2.innerHTML="请输入x₀的值"; return; } x0=parseFloat(Text1.value); if(isNaN(x0)){ spanr2.innerHTML="输入x₀的格式错误！"; return; } y0 = a + b * x0; yinc = ta * Math.sqrt(sig2 * (((n + 1.0) / n) + (x0 - xmean) * (x0 - xmean) / sigx2)); spanr2.innerHTML="y₀的预测区间是："+dtos(y0)+"±"+dtos(yinc)+",
"; spanr2.innerHTML += "也就是：(" + dtos(y0 - yinc) + ", " + dtos(y0 + yinc) + ")"; } function Button3_onclick() { var s; s="6,31,10,58,21,124,40,220,62,299,62,190,90,320,100,406,120,380"; TextArea1.value=s; spanresult.innerHTML=""; spanr2.innerHTML=""; divpre.style.display="none"; } function Button4_onclick() { var s; s="100,45,110,51,120,54,130,61,140,66,150,70,160,74,170,78,180,85,190,89"; TextArea1.value=s; spanresult.innerHTML=""; spanr2.innerHTML=""; divpre.style.display="none"; } function Radio1_onclick() { spanresult.innerHTML=""; spanr2.innerHTML=""; divpre.style.display="none"; } // JavaScript Document