cfs.py

# -*- coding: utf-8 -*-
"""CFS.ipynb

Automatically generated by Colaboratory.

Original file is located at
    https://colab.research.google.com/drive/1O4pGA5R6HkZKEFL_Q98v1pEZPxNUlhoJ
"""

import zipfile
import matplotlib.pyplot as plt
from math import log
import numpy as np
import pandas as pd
pd.options.mode.chained_assignment = None  # default='warn'

from scipy.stats import entropy
from sklearn.metrics import mutual_info_score

## PELT algorithm
#!pip install ruptures
import ruptures as rpt
import seaborn as sns
import collections


debug = False
events = []
workload = []
tag = 1
data = None

def change_point(data):
    ## offline change point detection: Pelt
    penalty = 2*np.log(len(data))*np.std(data)**2
    algo = rpt.Pelt(model = 'rbf').fit(data)
    bkps = algo.predict(pen = penalty)

    if debug:
      rpt.display(data,bkps,figsize=(10,6))
      plt.title('Change Point Detection: Pelt Search Method')
      plt.show()
      print(len(bkps))
      
    return bkps

# Load data from file
def load_data(fname):
  zf = zipfile.ZipFile(fname)

  global events
  events_data = pd.read_csv(zf.open("EVTS"))
  events = np.asarray(events_data['EVENTS'])
  ##events = np.asarray(events_data['EVENTS'])[1:12]
  print(events)

  raw_data = pd.read_csv(zf.open("DATA"))
  ##raw_data = pd.read_csv(zf.open("DATA")).iloc[0:100,[0,1,2,3,5,6,7,8,9,10,11,12,13,14,15]]
  print(raw_data)

  ##global workload for future usuage
  workload_data = pd.read_csv(zf.open('META'))  
  workload = np.asarray(workload_data['BENCHMARK'])

  if debug:
    print(raw_data)
    print(events)
    print(workload)
  
  return raw_data

##
def caclulate_expected_action(raw_data):
  ## add new columes 'PHASE' to save ave IPC during the phase
  ## Reward and Expected action
  raw_data.insert(raw_data.shape[1],'PHASE',0.0)
  raw_data.insert(raw_data.shape[1],'PHASE_ON',0.0)
  raw_data.insert(raw_data.shape[1],'PHASE_OFF',0.0)
  raw_data.insert(raw_data.shape[1],'REWARD',0.0)
  raw_data.insert(raw_data.shape[1],'EXPECTED',0)

  ## sort by CORE and then split it
  raw_data.sort_values('CORE',inplace=True)
  ## identify how many cores here
  core_counts = len(np.unique(raw_data['CORE']))
  split_data = np.vsplit(raw_data,core_counts)
  raw_data.sort_index(inplace=True)

  for i in range(core_counts):
    split_data[i].sort_values('ITER',inplace=True)
    data = np.asarray(split_data[i]['IPC'])
    mode = np.asarray(split_data[i]['MODE'])

    cps = change_point(data)
    if debug:
      print(split_data[i])
      print(cps)
    idx = 0
    for p in cps:
      #r = data[idx:p].mean()      
      r = 0.0
      val_on, cnt_on, r_on = 0.0,0,0.0
      val_off, cnt_off,r_off = 0.0,0,0.0
      for j in range(idx,p):
        if mode[j] == 1:
          cnt_on+=1
          val_on+=data[j]
        else:
          cnt_off+=1
          val_off+=data[j]
      if cnt_on:
        r_on = val_on / cnt_on
      else:
        r_on = 0.0
      if cnt_off:
        r_off = val_off / cnt_off
      else:
        r_off = 0.0
      if cnt_on + cnt_off:         
        r = (val_on + val_off) / (cnt_on + cnt_off)
      #print(r, r_on, r_off)
      while idx < p:
        raw_data['PHASE'][idx*core_counts+i] = r
        raw_data['PHASE_ON'][idx*core_counts+i] = r_on
        raw_data['PHASE_OFF'][idx*core_counts+i] = r_off
        idx += 1

  #print([raw_data['PHASE'],raw_data['PHASE_ON'],raw_data['PHASE_OFF']])
  ## Caculate reward and expected action

  i = 0
  while i  < (len(raw_data)-core_counts+1):
    reward = [0.0]*4
    rd = [0.0]*4
    rd_mode = [0.0]*4

    ## Caculate the reward 
    for j in range(core_counts):      
      if raw_data['PHASE'][i+j]:
        rd[j] = raw_data['IPC'][i+j]/raw_data['PHASE'][i+j]
      else:
        rd[j] = 0.0
      if raw_data['MODE'][i+j] == 1:
        if raw_data['PHASE_ON'][i+j]:
          rd_mode[j] =  raw_data['IPC'][i+j]/raw_data['PHASE_ON'][i+j]
        else:
          rd_mode[j] = 0.0
      else:
        if raw_data['PHASE_OFF'][i+j]:
          rd_mode[j] =  raw_data['IPC'][i+j]/raw_data['PHASE_OFF'][i+j]
        else:
          rd_mode[j] = 0.0

    val = 0.0
    for j in range(core_counts):
      val = 0.0
      for k in range(core_counts):
        if j == k:
          val += rd[k]
        else:
          val += rd_mode[k]
      
      reward[j] = val/core_counts -1

    for j in range(core_counts):
      raw_data['REWARD'][i+j] = reward[j]
      ## Caculate the expected action
      ## if reward is positive, expected action would be same as mode
      ## otherwise, oppsite of mode
      if (reward[j]) > 0 :
        raw_data['EXPECTED'][i+j] = raw_data['MODE'][i+j]
      else:
        raw_data['EXPECTED'][i+j] = 1- raw_data['MODE'][i+j]

    i += core_counts
 
  raw_data.sort_index(inplace=True)

  raw_data.drop(columns=['REWARD','IPC','MODE','PHASE','PHASE_ON','PHASE_OFF'],inplace=True)
  feature_counts = raw_data.shape[1] - 3

  if debug:
    print(raw_data)
    print(feature_counts)
    for i in range(feature_counts):
      plt.scatter(raw_data['FEATURE'+str(i)],raw_data['EXPECTED'],s=20)
      plt.title('Feature' + str(i) + ' : ' + events[i])
      plt.show()

    for i in range(feature_counts):
      #plt.scatter(raw_data['ITER'],raw_data['FEATURE'+str(i)])
      plt.scatter(raw_data.index,raw_data['FEATURE'+str(i)])
      plt.title('Feature' + str(i)+ ' : ' + events[i])
      plt.show()
    
  return raw_data

### Entropy Section ###
def get_probs(a):
    bases = collections.Counter([tmp_base for tmp_base in a])
    return [x/sum(bases.values()) for x in bases.values()]
        
def entropyd(a):
  return entropy(get_probs(a))

##mutual information
def conditional_entropy(f1, f2):
    return  entropyd(f1) - mutual_info_score(f1, f2) 

def information_gain(f1, f2):
    return entropyd(f1) - conditional_entropy(f1, f2)

def su_correlation(f1, f2):
    # calculate information gain of f1 and f2, t1 = ig(f1, f2)
    if (f1 == f2).all():
      return 1.0

    t1 = information_gain(f1, f2)
    t2 = entropyd(f1)
    t3 = entropyd(f2)
    su = 2.0 * t1 / (t2 + t3)

    return su

def gain(c,idx):
    '''    
    Return the information gain obtained by splitting a numeric attribute in two according to cut_point
    '''
    N = len(c)
    cl,cr = c[0:idx],c[(idx+1):N]
    return entropyd(c) - entropyd(cl)*idx/N - entropyd(cr)*(N-idx)/N

def get_gain(c):
  return [gain(c,i+1) for i in range(len(c)-2)]


def mdlpc_criterion(c,idx):
  N = len(c)
  k = len(np.unique(c))
  k1 = len(np.unique(c[0:idx+1]))
  k2 = len(np.unique(c[idx+1:N]))
  delta = log(3**k-2,2) - k*entropyd(c) + k1*(entropyd(c[0:idx+1]) + k2*entropyd(c[(idx+1):N]))
  gain_threshold = (log(N - 1, 2) + delta) / N
  return gain_threshold

def merits(select_sets,r_fc,r_ff):
    """
    This function calculates the merit of X given class labels y, where
    merits = (k * rcf)/sqrt(k+k*(k-1)*rff)
    rcf = (1/k)*sum(su(fi,y)) for all fi in X
    rff = (1/(k*(k-1)))*sum(su(fi,fj)) for all fi and fj in X
    """
    k = len(select_sets)
    if k == 0:
      return 0.0
      
    sum_r_fc,sum_r_ff = 0.0,0.0

    for s in select_sets:
      sum_r_fc += r_fc[s]

    cnt = 0
    for i in select_sets: 
      for j in select_sets:
        if i < j:
          cnt += 1
          sum_r_ff += r_ff[i][j]

    fc_mean = sum_r_fc/k
    if cnt !=0:
      ff_mean = sum_r_ff/(cnt)
    else:
      ff_mean = 0.0

    merit = (k * fc_mean) / np.sqrt(k + k*(k-1)*ff_mean)
    return merit

def forward_search(r_fc,r_ff):
  feature_counts = len(r_fc)

  feature_sets = [i for i in range(feature_counts)]
  forward_sets = []
  forward_dicts=[]
 
  preV = merits(forward_sets,r_fc,r_ff)

  while True:
    maxV,maxI = -float("inf"),-1
    
    for feature in range(feature_counts):
      ## cacaulate the maximum merits after adding a feature
      if feature not in forward_sets:
        forward_sets.append(feature)
        m=merits(forward_sets,r_fc,r_ff)
        if debug:
          print([m,maxV,preV])
        if (m > maxV):
          maxV = m
          maxI = feature
        forward_sets.remove(feature)
      
    if preV <= maxV and maxI !=-1:
        preV = maxV
        forward_dicts.append((events[maxI],maxV))
        forward_sets.append(maxI)
    else:
        break;

  print(forward_dicts)
  return forward_sets

def backward_search(r_fc,r_ff):
  feature_counts = len(r_fc)
  feature_sets = [i for i in range(feature_counts)]
  backward_sets = feature_sets.copy()
  backward_remove_order = []
  
  preV = merits(backward_sets,r_fc,r_ff)

  while True:
    maxV,maxI = -float("inf"),-1

    for feature in feature_sets:
      ## cacaulate the maximum merits after remving a feature 
        if feature in backward_sets:
          backward_sets.remove(feature)
          m=merits(backward_sets,r_fc,r_ff)
          if(m > maxV):
            maxV = m
            maxI = feature
          backward_sets.append(feature) 

    if preV <= maxV and maxI != -1:
      preV = maxV
      backward_remove_order.append((events[maxI],maxV))
      backward_sets.remove(maxI)
    else:
      break
  
  backward_rest_order = []
  for feature in backward_sets:
    tmp_sets = backward_sets.copy()
    tmp_sets.remove(feature)
    m=merits(tmp_sets,r_fc,r_ff)
    backward_rest_order.append([events[feature],m])

  print(backward_remove_order)
  print(backward_rest_order)
  return backward_sets

def partition(l,r):
  global data 

  tmpData = np.asarray(data['EXPECTED'][l:r]).copy()
  gains = get_gain(tmpData)
  if debug:
    plt.plot(gains)
    plt.show()
  maxI = np.argmax(gains)

  if debug:
    print([l,l+maxI+1,r])

  if (r-l)> 1 and gains[maxI] > mdlpc_criterion(tmpData, maxI) :
    partition(l,l+maxI+1)
    partition(l+maxI+1,r)
  else:
    if debug:
      print([l,r])
    for j in range(l,r):
      global tag
      data['G'][j]=tag
    tag += 1

## Final result
def main():

  raw_data = pd.DataFrame()

  for i in range(1): 
    fname = ""
    if i+5 < 10:
        fname = "./W00"+str(i+5)+".json.0.zip"
    else:
        fname = "./W0"+str(i+5)+".json.0.zip"

    tmp_raw_data = load_data(fname)  
    tmp_raw_data = caclulate_expected_action(tmp_raw_data)
    raw_data=raw_data.append(tmp_raw_data)
  
  n = len(raw_data)
  feature_counts = raw_data.shape[1] - 3
  g = [[0.0 for i in range(n)] for j in range(feature_counts)]
 
  for i in range(feature_counts):
    global tag
    tag = 1
    feature_name = 'FEATURE' + str(i)
    if debug:
      print(feature_name+ ' : ' + events[i])

    raw_data.sort_values(feature_name,inplace=True)
    global data
    data = pd.DataFrame({'EXPECTED':raw_data['EXPECTED'],'G':-1,'Iter':raw_data.index})
    
    data.index = range(n)

    partition(0,n)

    data.sort_values('Iter',inplace=True)
    g[i] = np.asarray(data['G']).copy()
    data.drop(data.index, inplace=True)

  if debug:
    for i in range(feature_counts):
      print(g[i])
      plt.ylim(ymax=4,ymin=0)
      plt.scatter(raw_data['FEATURE'+str(i+feature_index)],g[i])
      plt.title('Feature' + str(i) + ' : ' + events[i])
      plt.show()

  raw_data.sort_index(inplace=True)
  c = np.asarray(raw_data['EXPECTED']).copy()

  r_fc = [0.0 for i in range(feature_counts)]
  r_ff = [[0.0 for i in range(feature_counts)] for j in range(feature_counts)]

  for i in range(feature_counts):
    r_fc[i]=su_correlation(g[i],c)
    for j in range(feature_counts):
        r_ff[i][j] = su_correlation(g[i],g[j])

  label = events
  pd_r_ff = pd.DataFrame(r_ff,columns=label,index=label)
  pd_r_ff.to_csv("r_ff.csv", index=True)
  sns.heatmap(pd_r_ff,cmap="RdBu_r")

  plt.savefig('feature.png', dpi=300)
  plt.close()

  if debug:
    print(r_fc)
    for i in range(feature_counts):
      print(r_ff[i]) 

  fsets = forward_search(r_fc,r_ff) 
  bsets = backward_search(r_fc,r_ff)
  res = np.union1d(fsets,bsets)
  return [(i,events[int(i)]) for i in res]

print(main())