In [1]:
import numpy as np
import pandas
import scipy, scipy.spatial
import sklearn
import sys

from matplotlib import pyplot as plt

%matplotlib inline

In [2]:
y = pandas.read_table("~/Downloads/data/ml/label_train.txt", sep=" ", dtype='int', header=None)

ndim= 900
y.head()


Out[2]:
0
0 161
1 163
2 56
3 119
4 138

In [3]:
ymin = 157
ysplit = 159
ymax = 159

In [4]:
np.unique(y[0], return_counts=True)


Out[4]:
(array([  1,   2,   3,   4,   5,   6,   7,   8,   9,  10,  11,  12,  13,
         14,  15,  16,  17,  18,  19,  20,  21,  22,  23,  24,  25,  26,
         27,  28,  29,  30,  31,  32,  33,  34,  35,  36,  37,  38,  39,
         40,  41,  42,  43,  44,  45,  46,  47,  48,  49,  50,  51,  52,
         53,  54,  55,  56,  57,  58,  59,  60,  61,  62,  63,  64,  65,
         66,  67,  68,  69,  70,  71,  72,  73,  74,  75,  76,  77,  78,
         79,  80,  81,  82,  83,  84,  85,  86,  87,  88,  89,  90,  91,
         92,  93,  94,  95,  96,  97,  98,  99, 100, 101, 102, 103, 104,
        105, 106, 107, 108, 109, 110, 111, 112, 113, 114, 115, 116, 117,
        118, 119, 120, 121, 122, 123, 124, 125, 126, 127, 128, 129, 130,
        131, 132, 133, 134, 135, 136, 137, 138, 139, 140, 141, 142, 143,
        144, 145, 146, 147, 148, 149, 150, 151, 152, 153, 154, 155, 156,
        157, 158, 159, 160, 161, 162, 163, 164]),
 array([  1263,   1261,   1255,   1256,   1252,   1235,   1240,   1264,
          1256,   1281,   1245,   1278,   1278,   1253,   1255,   1255,
          1291,   1277,   1308,   1285,   1322,   1309,   1318,   1322,
          1327,   1339,   1361,   1361,   1335,   1396,   1359,   1393,
          1373,   1356,   1398,   1416,   1386,   1398,   1396,   1404,
          1430,   1398,   1416,   1406,   1420,   1445,   1433,   1445,
          1454,   1451,   1481,   1482,   1477,   1474,   1478,   1486,
          1512,   1492,   1557,   1557,   1548,   1530,   1574,   1582,
          1606,   1611,   1666,   1650,   1704,   1739,   1735,   1743,
          1728,   1796,   1737,   1810,   1822,   1864,   1847,   1838,
          1857,   1913,   1910,   1917,   2006,   1992,   2033,   2063,
          2072,   2063,   2096,   2128,   2134,   2206,   2215,   2212,
          2258,   2279,   2287,   2319,   2356,   2435,   2438,   2491,
          2486,   2485,   2502,   2555,   2594,   2629,   2575,   2587,
          2777,   2875,   2897,   2884,   2978,   3087,   3179,   3368,
          3388,   3421,   3409,   3453,   3536,   3586,   3615,   3696,
          3821,   3802,   3934,   4059,   4069,   4253,   4819,   4939,
          5038,   5259,   5310,   6080,   6487,   6623,   7256,   8279,
          9069,   9221,   9707,   9998,  10557,  10645,  11484,  12382,
         12858,  16548,  18562,  21943,  30679,  34092,  45439,  60513,
         64478,  65211,  92241, 130122]))

In [5]:
yuniq,ycount = np.unique(y[0], return_counts=True)

print(np.sum(ycount[np.where(np.in1d(yuniq, range(ymin, ysplit)))[0]]))
print(np.sum(ycount[np.where(np.in1d(yuniq, range(ysplit, ymax+1)))[0]]))


64771
45439

In [6]:
import pickle

cstat = pickle.load(open( "../data/sum_features.dat", "rb" ) )

In [7]:
### Calclulate Standardized Mean Difference Between Classes

def calStandMeanDiff(y, cstat, yneg, ypos):
    sx  = np.zeros(shape=ndim, dtype=float)
    ssx = np.zeros(shape=ndim, dtype=float)


    n1 = np.sum(np.in1d(y, yneg))
    n2 = np.sum(np.in1d(y, ypos))
    sys.stderr.write("Number of samples in NegClass: %d and PosClass: %d \n"%(n1, n2))

    for yi in yneg:
        sx += cstat[yi][0]
        ssx += cstat[yi][1]
    r1_mean = sx / float(n1)
    r1_var = (ssx - 2*sx*r1_mean + r1_mean**2) / float(n1)

    tot_mean = sx
    tot_var  = ssx
    
    sx  = np.zeros(shape=ndim, dtype=float)
    ssx = np.zeros(shape=ndim, dtype=float)
    for yi in ypos:
        sx += cstat[yi][0]
        ssx += cstat[yi][1]
    r2_mean = sx / float(n2)
    r2_var = (ssx - 2*sx*r2_mean + r2_mean**2) / float(n2)

    tot_mean += sx
    tot_var  += ssx
    tot_mean = tot_mean / float(n1 + n2)
    tot_var  = (tot_var - 2*tot_var*tot_mean + tot_mean**2) / float(n1 + n2)

    rdiff = (r1_mean - r2_mean) / np.sqrt(tot_var)

    return (rdiff)


## unit test:
mean_test = calStandMeanDiff(y, cstat, np.arange(ymin,ysplit), np.arange(ysplit, ymax+1)) 
print(np.sum(mean_test > 0.001))


175
Number of samples in NegClass: 64771 and PosClass: 45439 

Classify items belonging to first half (1) Second half (-1)

Finding Good Features


In [8]:
rdiff = calStandMeanDiff(y, cstat, np.arange(ymin,ysplit), np.arange(ysplit, ymax+1))


## Good Features:
goodfeatures = np.where(rdiff > 0.001)[0]

goodfeatures


Number of samples in NegClass: 64771 and PosClass: 45439 
Out[8]:
array([  6,  19,  35,  39,  40,  59,  73,  75,  77,  80,  87,  91,  94,
       108, 109, 114, 125, 137, 138, 140, 151, 153, 155, 159, 166, 175,
       178, 180, 185, 187, 190, 193, 197, 202, 210, 214, 218, 236, 239,
       243, 245, 247, 249, 257, 259, 263, 265, 267, 272, 273, 277, 278,
       281, 282, 286, 295, 298, 307, 308, 328, 329, 332, 333, 334, 340,
       341, 344, 347, 348, 350, 352, 353, 355, 357, 359, 366, 373, 375,
       387, 391, 399, 407, 412, 416, 418, 422, 430, 434, 439, 443, 444,
       448, 457, 460, 461, 462, 467, 489, 491, 500, 501, 506, 517, 526,
       527, 530, 531, 538, 541, 542, 544, 556, 561, 564, 566, 571, 573,
       579, 581, 584, 594, 600, 608, 613, 630, 652, 657, 665, 671, 676,
       680, 686, 696, 705, 710, 712, 714, 716, 720, 722, 736, 738, 746,
       747, 749, 760, 766, 770, 772, 773, 774, 779, 794, 798, 799, 802,
       805, 812, 815, 816, 817, 825, 829, 853, 855, 859, 863, 864, 869,
       873, 877, 881, 886, 890, 892])

Read a Random Sample


In [9]:
def readRandomSample(data_fname, y, size, goodfeat=None, acc_miny=None, acc_maxy=None):
    """ Read a random sample
    """
    if goodfeat is None:
        goodfeat = np.arange(ndim)
    Xsub = np.empty(shape=(size,goodfeat.shape[0]), dtype=float)
    ysub = np.zeros(shape=size, dtype=int)

    if acc_miny is None:
        acc_miny = np.min(y)
    if acc_maxy is None:
        acc_maxy = np.max(y)
        
    #yuniq, ycount = np.unique(y, return_counts=True)
    #tot_acceptable = np.sum(ycount[np.where((yuniq >= acc_miny) & (yuniq <= acc_maxy))[0]])
    
    acceptable_indx = np.where((y>=acc_miny) & (y<=acc_maxy))[0]
    assert(acceptable_indx.shape[0] > size)
    choice_indx = np.sort(np.random.choice(acceptable_indx, size, replace=False))
    #print(choice_indx.shape)
    #sys.stderr.write("Total Accetables: --> %d"%(tot_acceptable))
    
    #proba = 1.0 - size/float(tot_acceptable)
    
        
    with open(data_fname, 'r') as fp:
        n = 0
        nf = 0
        for line in fp:
#            if (y[n] >= acc_miny and y[n]<=acc_maxy):
#                if np.random.uniform(low=0, high=1) > proba and nf < size:
            if nf < size:
                if n == choice_indx[nf]:
                    line = line.strip().split()
                    ix = -1
                    for i,v in enumerate(line):
                        if np.any(goodfeat == i):
                            ix += 1
                            Xsub[nf,ix] = int(v)
                    ysub[nf] = y[n]

                    nf += 1
            n += 1
    return(Xsub, ysub)

In [10]:
## unit testing readRandomSample()
gf_test = goodfeatures
Xsub, ysub = readRandomSample('/home/vahid/Downloads/data/ml/data_train.txt', y[0], \
                              size=2000, goodfeat=gf_test, acc_miny=ymin, acc_maxy=ymax)

print(Xsub.shape)
print(np.unique(ysub))


(2000, 175)
[157 158 159]

In [11]:
### Performance Evaluation
def evalPerformance(ytrue, ypred):
    tp = np.sum(ypred[np.where(ytrue ==  1)[0]] == 1)
    fp = np.sum(ypred[np.where(ytrue == -1)[0]] == 1)
    tn = np.sum(ypred[np.where(ytrue == -1)[0]] == -1)
    fn = ytrue.shape[0]-(tp+fp+tn)
    #sys.stderr.write('%d %d %d %d\n'%(tp,fp,tn,fn))
    prec = tp / float(tp + fp)
    recall  = tp / float(tp + fn)
    f1score = 2*tp/float(2*tp + fp + fn)

    return (prec, recall, f1score)

In [12]:
Xsub, ysub = readRandomSample('/home/vahid/Downloads/data/ml/data_train.txt', y[0], size=20000, \
                              goodfeat=goodfeatures, acc_miny=ymin, acc_maxy=ymax)

assert(np.sum(ysub < ymin) == 0)
assert(np.sum(ysub > ymax) == 0)
ysub[np.where(ysub < ysplit)[0]] = -1
ysub[np.where(ysub >= ysplit)[0]] =  1

print(np.sum(ysub == -1), np.sum(ysub==1))

#Xsub = Xsub[:, goodfeatures]
x_std = np.std(Xsub, axis=0)
print(np.where(x_std < 0.0001))
Xsub = (Xsub - np.mean(Xsub, axis=0)) / np.std(Xsub, axis=0)

Xsub.shape


(11785, 8215)
(array([], dtype=int64),)
Out[12]:
(20000, 175)

Grid-Search (coarse)


In [13]:
import sklearn.svm

ntot = Xsub.shape[0]
tr_idx = np.random.choice(ntot, size=ntot/2, replace=False)
ts_idx = np.setdiff1d(np.arange(ntot), tr_idx, assume_unique=True)
yts = ysub[ts_idx]

for c in [0.01, 0.1, 1.0, 5.0, 10.0]:
    for gm in [0.001, 0.01, 0.1, 1.0, 5.0]:
        clf = sklearn.svm.SVC(C=c, kernel='rbf', gamma=gm)
        clf.fit(Xsub[tr_idx, :], ysub[tr_idx])
        ypred = clf.predict(Xsub[ts_idx, :])
        prec, recall, f1score = evalPerformance(yts, ypred)
        print ("C=%.4f Gamma=%.4f  ==> Prec:%.3f  Recall:%.3f  F1Score:%.3f"%(c, gm, prec, recall, f1score))


C=0.0100 Gamma=0.0010  ==> Prec:nan  Recall:0.000  F1Score:0.000
C=0.0100 Gamma=0.0100  ==> Prec:0.711  Recall:0.428  F1Score:0.534
C=0.0100 Gamma=0.1000  ==> Prec:0.831  Recall:0.030  F1Score:0.058
C=0.0100 Gamma=1.0000  ==> Prec:nan  Recall:0.000  F1Score:0.000
C=0.0100 Gamma=5.0000  ==> Prec:nan  Recall:0.000  F1Score:0.000
C=0.1000 Gamma=0.0010  ==> Prec:0.666  Recall:0.644  F1Score:0.655
C=0.1000 Gamma=0.0100  ==> Prec:0.724  Recall:0.534  F1Score:0.614
C=0.1000 Gamma=0.1000  ==> Prec:0.806  Recall:0.124  F1Score:0.215
C=0.1000 Gamma=1.0000  ==> Prec:0.778  Recall:0.002  F1Score:0.003
C=0.1000 Gamma=5.0000  ==> Prec:nan  Recall:0.000  F1Score:0.000
C=1.0000 Gamma=0.0010  ==> Prec:0.681  Recall:0.671  F1Score:0.676
C=1.0000 Gamma=0.0100  ==> Prec:0.726  Recall:0.616  F1Score:0.666
C=1.0000 Gamma=0.1000  ==> Prec:0.784  Recall:0.225  F1Score:0.349
C=1.0000 Gamma=1.0000  ==> Prec:0.796  Recall:0.010  F1Score:0.019
C=1.0000 Gamma=5.0000  ==> Prec:0.600  Recall:0.001  F1Score:0.003
C=5.0000 Gamma=0.0010  ==> Prec:0.688  Recall:0.690  F1Score:0.689
C=5.0000 Gamma=0.0100  ==> Prec:0.712  Recall:0.641  F1Score:0.675
C=5.0000 Gamma=0.1000  ==> Prec:0.787  Recall:0.232  F1Score:0.358
C=5.0000 Gamma=1.0000  ==> Prec:0.836  Recall:0.012  F1Score:0.025
C=5.0000 Gamma=5.0000  ==> Prec:0.636  Recall:0.002  F1Score:0.003
C=10.0000 Gamma=0.0010  ==> Prec:0.690  Recall:0.698  F1Score:0.694
C=10.0000 Gamma=0.0100  ==> Prec:0.704  Recall:0.629  F1Score:0.664
C=10.0000 Gamma=0.1000  ==> Prec:0.784  Recall:0.230  F1Score:0.356
C=10.0000 Gamma=1.0000  ==> Prec:0.836  Recall:0.012  F1Score:0.025
C=10.0000 Gamma=5.0000  ==> Prec:0.636  Recall:0.002  F1Score:0.003

In [14]:
import sklearn.svm

ntot = Xsub.shape[0]
tr_idx = np.random.choice(ntot, size=ntot/2, replace=False)
ts_idx = np.setdiff1d(np.arange(ntot), tr_idx, assume_unique=True)
yts = ysub[ts_idx]

for c in [8, 10, 12, 15, 20]:
    for gm in [0.0005, 0.0008, 0.001, 0.0015, 0.005, 0.008]:
        clf = sklearn.svm.SVC(C=c, kernel='rbf', gamma=gm)
        clf.fit(Xsub[tr_idx, :], ysub[tr_idx])
        ypred = clf.predict(Xsub[ts_idx, :])
        prec, recall, f1score = evalPerformance(yts, ypred)
        print ("C=%.4f Gamma=%.4f  ==> Prec:%.3f  Recall:%.3f  F1Score:%.3f"%(c, gm, prec, recall, f1score))


C=8.0000 Gamma=0.0005  ==> Prec:0.662  Recall:0.702  F1Score:0.681
C=8.0000 Gamma=0.0008  ==> Prec:0.669  Recall:0.707  F1Score:0.687
C=8.0000 Gamma=0.0010  ==> Prec:0.671  Recall:0.710  F1Score:0.690
C=8.0000 Gamma=0.0015  ==> Prec:0.674  Recall:0.724  F1Score:0.698
C=8.0000 Gamma=0.0050  ==> Prec:0.661  Recall:0.703  F1Score:0.681
C=8.0000 Gamma=0.0080  ==> Prec:0.675  Recall:0.672  F1Score:0.674
C=10.0000 Gamma=0.0005  ==> Prec:0.663  Recall:0.706  F1Score:0.684
C=10.0000 Gamma=0.0008  ==> Prec:0.669  Recall:0.709  F1Score:0.689
C=10.0000 Gamma=0.0010  ==> Prec:0.672  Recall:0.717  F1Score:0.693
C=10.0000 Gamma=0.0015  ==> Prec:0.673  Recall:0.729  F1Score:0.700
C=10.0000 Gamma=0.0050  ==> Prec:0.656  Recall:0.702  F1Score:0.678
C=10.0000 Gamma=0.0080  ==> Prec:0.672  Recall:0.666  F1Score:0.669
C=12.0000 Gamma=0.0005  ==> Prec:0.665  Recall:0.707  F1Score:0.685
C=12.0000 Gamma=0.0008  ==> Prec:0.670  Recall:0.712  F1Score:0.690
C=12.0000 Gamma=0.0010  ==> Prec:0.671  Recall:0.719  F1Score:0.694
C=12.0000 Gamma=0.0015  ==> Prec:0.671  Recall:0.733  F1Score:0.701
C=12.0000 Gamma=0.0050  ==> Prec:0.653  Recall:0.701  F1Score:0.676
C=12.0000 Gamma=0.0080  ==> Prec:0.669  Recall:0.661  F1Score:0.665
C=15.0000 Gamma=0.0005  ==> Prec:0.666  Recall:0.711  F1Score:0.688
C=15.0000 Gamma=0.0008  ==> Prec:0.671  Recall:0.719  F1Score:0.694
C=15.0000 Gamma=0.0010  ==> Prec:0.672  Recall:0.724  F1Score:0.697
C=15.0000 Gamma=0.0015  ==> Prec:0.670  Recall:0.734  F1Score:0.700
C=15.0000 Gamma=0.0050  ==> Prec:0.650  Recall:0.696  F1Score:0.672
C=15.0000 Gamma=0.0080  ==> Prec:0.663  Recall:0.657  F1Score:0.660
C=20.0000 Gamma=0.0005  ==> Prec:0.667  Recall:0.710  F1Score:0.688
C=20.0000 Gamma=0.0008  ==> Prec:0.670  Recall:0.721  F1Score:0.695
C=20.0000 Gamma=0.0010  ==> Prec:0.670  Recall:0.730  F1Score:0.699
C=20.0000 Gamma=0.0015  ==> Prec:0.664  Recall:0.735  F1Score:0.698
C=20.0000 Gamma=0.0050  ==> Prec:0.646  Recall:0.693  F1Score:0.669
C=20.0000 Gamma=0.0080  ==> Prec:0.659  Recall:0.648  F1Score:0.653

In [ ]: