C++でパーセプトロン

練習としてパーセプトロン（線形学習器）をC++で実装した。
動くには動くが多分バグがあると思う。検証しろよ。

アルゴリズムはCristianini, Shawe-Taylor著 大北訳の「サポートベクターマシン入門」（赤本）を参考にした。
データフォーマットはLIBSVMのものが使えるようにしました。

入力は訓練セットとテストセット、出力はテストセットの正解率など。
なお線形分離不可能なデータセットは学習できません（無限ループになる）。

多分まともにプログラムできる人からしたらゴミみたいなコードだと思います
もっとシンプルに美しく組めるんだろうなぁ。

できれば2次元の場合の可視化をしたい。

#include <cstdlib>
#include <iostream>
#include <numeric>
#include <vector>
#include <math.h>

using namespace std;

struct str{
    int label;
    double b;
    vector<int> index;
    vector<double> value;
};

double p=1;    //learning rate
vector <str> trainset;
vector <str> testset;
str classfier;


void input(char* fname, vector<str>* instance){
    int n=0,c=0;
    double v;    
    str st;
    FILE *fin=fopen(fname, "r");
    
    while(1){
        int i=0;
        fscanf(fin,"%d",&st.label);
        while(1){
            do{
                c=fgetc(fin);
            if(c=='\n') goto out;
            }while(isspace(c));
            ungetc(c,fin);
            fscanf(fin,"%d:%lf",&n,&v);
            st.index.push_back(n);
            st.value.push_back(v);
            ++i;
        }
        out:        
        instance->push_back(st);
        for(int j=0;j<i;++j){
            st.index.pop_back();
            st.value.pop_back();
        }        
        c=fgetc(fin);
        if(c==EOF) break;
        ungetc(c,fin);
    }
}

double inner_product(str* x, str* y){
    int i,j;
    double sum=0;
    for(i=0;i<x->index.size();++i){
        for(j=i; j<y->index.size(); ++j){
            if(x->index[i]==y->index[j]){
                sum += x->value[i] * y->value[j];
            }            
        }
    }
    return sum;
}

void vector_sum(str* x, str* y){
    int i=0,j=0;
    str tmp_str;
    
    while(1){
        if(x->index[i]==y->index[j]){
            tmp_str.index.push_back(x->index[i]);
            tmp_str.value.push_back(p*x->label*x->value[i]+y->value[j]); 
            ++i; ++j;
        }
        else if(x->index[i] < y->index[j]){
            tmp_str.index.push_back(x->index[i]);
            tmp_str.value.push_back(p*x->label*x->value[i]);
            ++i;
        }
        else{
            tmp_str.index.push_back(y->index[j]);
            tmp_str.value.push_back(y->value[j]);
            ++j;
        }
        
        if(i>x->index.size()-1 || j>y->index.size()-1){
            if(j<y->index.size()){
                while(j<y->index.size()){
                    tmp_str.index.push_back(y->index[j]);
                    tmp_str.value.push_back(y->value[j]);
                    ++j;
                }
            }
            else if(i<x->index.size()){
                while(i<x->index.size()){
                    tmp_str.index.push_back(x->index[i]);
                    tmp_str.value.push_back(p*x->label*x->value[i]);
                    ++i;
                }
            }
            break;
        }
    }    
    classfier = tmp_str;
}

double function(str *w, str *x){
    return (x->label * (inner_product(w,x) + w->b));
}

void perceptrons(){
    int i,update=0,it=0;
    double R=0, tmp;
    classfier.b=0;
    
    for(i=0;i<trainset.size();++i){
        tmp=inner_product(&trainset[i],&trainset[i]);
        if(tmp>R) R=tmp;
    }
    R=sqrt(R);
    cout << "R=" << R << endl;
    classfier.index.push_back(trainset[0].index[0]);
    classfier.value.push_back(0);
    
    do{
        update=0;
        cout << "it : " << it <<endl;
        for(i=0;i<trainset.size();++i){
            cout << "f=" << function(&classfier, &trainset[i]) << " ";
            cout << "b=" << classfier.b << endl; 
            if(function(&classfier, &trainset[i]) <= 0){
                vector_sum(&trainset[i], &classfier);
                classfier.b=classfier.b+p*trainset[i].label*R*R;
                ++update;
            }
        } 
        ++it;
    }while(update!=0);
    
}

void test(){
    double acc=0;
    for(int i=0;i<testset.size();++i){
        if(function(&classfier,&testset[i]) > 0){
            ++acc;
        }
    }
    cout << "accuraccy : " << acc / testset.size() << endl;
}

int main(int argc, char** argv) {
    if(argc!=2){
        cout << "input trainfile testfile" <<endl;
    }
    input(argv[1],&trainset);
    input(argv[2],&testset);
    
    cout << "learn " << argv[1] << endl;
    
    perceptrons();
    
    cout << endl;
    cout << "b=" << classfier.b << endl;
    for(int i=0;i<classfier.index.size();++i){
        cout << classfier.index[i] << ":" << classfier.value[i] << endl;
    }
    
    cout << endl;
    cout<< "test " << argv[2] << endl;
    test();

    return 0;
}