# install packages

from matplotlib import pyplot as plt
import numpy as np
from qpsolvers import solve_qp
from matplotlib.patches import Arc

# assign input data

left_points = [2, 3]
right_points = [6, 6.5, 8.3]
addtl_points = [10.5]

# define a customized phi transformation

offset = np.average(7) 
def transform(x):
	return (x - offset)**2

# define the kernel function for the phi transformation

def decomp_k(xm, xn):
	return xm*xn + transform(xm)*transform(xn)

# alternate kernel function (linear)

def linear_k(xm, xn):
	return xm*xn

# alternate kernal function (polynomial)

def poly_k(xm, xn):
	c = 1
	d = 2	
	return (xm*xn + c)**d

# alternate kernel function (exponential)

def rbf_k(xm, xn):
	sigma = 1
	return np.exp(-(xm - xn)**2/(2*sigma**2))

# alternate kernel function (minimum)

def min_k(xm,xn):
	return min(xm,xn)

# select the preferred kernel

def k(xm, xn):
	return decomp_k(xm, xn)

# assemble the input data and classifications

x = left_points + right_points + addtl_points
y = np.concatenate((np.full((len(left_points), 1), -1),
	np.full((len(right_points), 1), 1),
	np.full((len(addtl_points), 1), -1)))

# assemble the matrices used for quadratic programming

P = np.zeros((len(x), len(x)))
for i in range(len(x)):
	for j in range(len(x)):
		P[i,j] = y[i]*y[j]*k(x[i], x[j])
print(P)
q = np.asarray([-1.]*len(x)).T
G = -1*np.identity(len(x))
h = np.asarray([0.]*len(x))		
A = np.array(y).T
b = np.array([0])		

P = P+1e-3*np.eye(len(x)) # stabilizes solution

# print results

solution = solve_qp(P, q, G, h, A, b) 
print('Solution: ', solution)
support_vectors = [i+1 for i,x in enumerate(G@np.asarray(solution)		
) if abs(x) > 1e-6]
print('Support vector indices: ', support_vectors)
def kprime(x0):
	return np.diag([k(x0, i) for i in x])
b = np.average([y[i - 1]-solution.T@kprime(x[i - 1])@y for i in support_vectors])
print('Bias: ', b)
for idx,i in enumerate(x):
	print(y[idx], solution.T@kprime(i)@y + b)

# find division points by checking where the classification output changes

def classification(x):
	return np.sign(solution.T@kprime(x)@y + b)[0]
class_input_array = np.linspace(0, 12, 5000)
class_output = [classification(x) for x in class_input_array]
zero_crossings = np.where(np.diff(np.sign(class_output)))[0]
(div_1, div_2) = class_input_array[zero_crossings]
print('Divisions: ', div_1, div_2)	

# plot the results

font = {'size':16}
tickfont = {'size':10}
fig, ax = plt.subplots()
plt.rcParams['axes.facecolor'] = 'g'
plt.rcParams['savefig.facecolor'] = '#c0c0c0'
plt.rc('font', **font)
marker_style = dict(linestyle='none', markersize = 10, markeredgecolor = 'k')

ax.set_aspect(1)
ax.plot(left_points,np.zeros(len(left_points)), marker = "^",
	markerfacecolor = 'red', **marker_style)
ax.plot(right_points, np.zeros(len(right_points)), marker = "s",
	markerfacecolor = 'blue', **marker_style)
ax.plot(addtl_points,np.zeros(len(addtl_points)), marker = "^",
	markerfacecolor = 'red', **marker_style)

for i in [div_1,div_2]:
	ax.annotate('', xy = (i, -0.6), xytext = (i, 0.6),
	arrowprops = {'arrowstyle':'-', 'lw':5, 'color':'0'},
	va = 'center')
arc = Arc(((div_1 + div_2)/2, 0.6), div_2 - div_1, 3, linestyle = '--',
	angle = 0.0, theta1 = 0, theta2 = 180)		
ax.add_artist(arc)
plt.gca().axes.get_yaxis().set_visible(False)
ax.annotate('', xy = (0,1.9), xytext = (0, -1.9),
arrowprops = {'arrowstyle':'-', 'lw':1, 'facecolor':'k'},
va = 'center')
fig.set_facecolor('#c0c0c0')
plt.xlim([-1, 13])
plt.ylim([-2.5, 2.5])
plt.xticks(np.arange(2,11,2), **tickfont)
ax.annotate('', xy = (12, 0), xytext = (-.05, 0),
arrowprops = {'arrowstyle':'->', 'lw':1, 'facecolor':'k'},
va = 'center', zorder = -1)
ax.annotate('$x$', xy = (0, 0), xytext = (12.2, 0), va = 'center') 
ax.text(1.4, 1, '-1 class', fontsize = '14')
ax.annotate('+1 class', xy = (0,0), xytext = ((div_1 + div_2)/2, 0.8),
	ha = 'center', fontsize = '14')
plt.tick_params(axis='x', which = 'major', labelsize = 10)
for i in support_vectors:
	ax.plot(x[i - 1], 0, marker= "o", linestyle = 'none',
	markersize = 25, markerfacecolor = 'None', markeredgecolor = 'k')

plt.box(on = None)
ax.spines['bottom'].set_position('zero')
ax.spines['bottom'].set_color('none')
ax.spines['left'].set_position('zero')
ax.spines['left'].set_color('none')
ax.spines['left'].set_linewidth(1)
ax.spines['bottom'].set_linewidth(1)

fig.savefig('graph.png', bbox_inches = 'tight', pad_inches = -0.1)
plt.show()