分散と相関係数(2) - rscのブログ

　分散と相関係数などを求めるプログラムをCで書いてみました。ついでに、前回の自分の質問で分かったコンソールマクロも付けてみました。

● correl.c

/*
 *  分散と相関係数 correl.c
 */

#include <stdio.h>
#include <conio.h>
#include <math.h>

/* Borland C++ 5.5.1 for Win32用 <conio.h> が必要 */
#define CLS()           clrscr()
#define LOCATE(x,y)     gotoxy((x)+1,(y)+1)         /* (0-79,0-24) */
#define COLOR(c)        textcolor(((c)%8) ? ((c)%8)+8 : 0)  /* cprintf()用 */
     /* 0:黒; 1:青; 2:緑; 3:水; 4:赤; 5:紫; 6:黄; 7:白 */
#define POS_0()         (wherex()-1)
#define CSRLIN()        (wherey()-1)

/* BASICﾗｲｸな倍精度用出力ﾏｸﾛ ; p: 先行 ｽﾄﾘﾝｸﾞ; q: 後続 ｽﾄﾘﾝｸﾞ */
#define PRINT_D(p,x,q)  printf("%s% .16g%s",p,(x),q)
#define ARRAY_SIZE(a)   (sizeof(a) / sizeof(a[0]))  /* 配列の要素数を求めるマクロ */

/* 関数名はExcelの関数に由来 */
double average_x(int, double *);                    /* 平均値 */
double average_xy(int, double *, double *);         /* 平均値 xy */
double varp(int, double *);                         /* 分散 母数 n */
double stdevp(int, double *);                       /* 標準偏差 母数 n */
double correl(int, double *,double *);              /* 相関係数 */

int main(void)
{
    double x[]={2.0,4.0,5.0};
    double y[]={3.0,3.0,4.0};
    int n=ARRAY_SIZE(x);

    CLS();
    PRINT_D("xの平均値    ＝",average_x(n,x),   "\n");
    PRINT_D("yの平均値    ＝",average_x(n,y),   "\n");
    PRINT_D("xyの平均値   ＝",average_xy(n,x,y),"\n");
    PRINT_D("xの分散      ＝",varp(n,x),        "\n");
    PRINT_D("yの分散      ＝",varp(n,y),        "\n");
    PRINT_D("xの標準偏差  ＝",stdevp(n,x),      "\n");
    PRINT_D("yの標準偏差  ＝",stdevp(n,y),      "\n");
    PRINT_D("x,yの相関係数＝",correl(n,x,y),    "\n");
    
    return(0);
}

double average_x(int n, double *x)
{
    int i=n;
    double s=0.0;
    
    while(i--)  /* n回ループ */
        s += *x++;
    if(n>0) return(s/(double)n);
    else    return(0.0);
}

double average_xy(int n, double *x, double *y)
{
    int i=n;
    double s=0.0;
    
    while(i--)
        s += *x++ * *y++;
    if(n>0) return(s/(double)n);
    else    return(0.0);
}

double varp(int n, double *x)
{
    int i=n;
    double s=0.0;
    double m=average_x(n,x);
    
    while(i--)
        s += *x * *x++;
    if(n>0) return(s/(double)n-m*m);
    else    return(0.0);

}

double stdevp(int n, double *x)
{
    return(sqrt(varp(n,x)));
}

double correl(int n, double *x,double *y)
{
    double mx=average_x(n,x), my=average_x(n,y);
    double mxy=average_xy(n,x,y);
    double sx=stdevp(n,x), sy=stdevp(n,y);
    
    if(sx*sy!=0.0) return((mxy-mx*my)/(sx*sy));
    else           return(0.0); 
}

●実行結果

xの平均値    ＝ 3.666666666666667
yの平均値    ＝ 3.333333333333333
xyの平均値   ＝ 12.66666666666667
xの分散　    ＝ 1.555555555555557
yの分散  　  ＝ 0.2222222222222212
xの標準偏差  ＝ 1.247219128924647
yの標準偏差  ＝ 0.4714045207910306
x,yの相関係数＝ 0.7559289460184548

※参考URL
http://q.hatena.ne.jp/1278083887
●分散と相関係数