Введение в машинное обучение.
В данной заметке рассматривются такие понятия как классы и иные библиотек относящихся к анализу данных.
Это предварительная версия! Любые замечания приветсвуются.
%matplotlib inline
import numpy as np
import matplotlib.pyplot as plt
import scipy.stats as models
a = np.random.normal(size=100)
a_mean = np.mean( a )
a_std = np.std( a )
a_mean, a_std
b = np.random.normal(size=100)
b_mean = np.mean( b )
b_std = np.std( b )
b_mean, b_std
class stats:
def calc( self, x ):
self.mean = np.mean( x )
self.std = np.std( x )
first = stats()
stats.calc( first, a )
print( first.mean, first.std ) # Обе статистики теперь вместе.
second = stats()
stats.calc( second, b )
print( second.mean, second.std ) # Обе статистики теперь вместе.
Классы на то и классы, что вызов можно сократить:
first.calc( a )
second.calc( b )
Печать тоже можно сделать методом
class stats:
def calc( self, x ):
self.mean = np.mean( x )
self.std = np.std( x )
def pri( self ):
print( self.mean, self.std)
first.pri()
first = stats()
stats.calc( first, a )
first.pri()
Добавление элемента в объект классв.
first.qq = 3 # Добавили переменную в класс.
first.qq
second.qq
third = stats()
third.pri() # Отсутвуют члены класса, так как они не были ещё созданы.
stats.__dict__
first.__dict__
second.__dict__
from types import MethodType
def get_mean( self ):
return self.mean
first.getmm=get_mean
first.getmm()
first.__dict__
first.getmm(first)
first.getmm = MethodType(get_mean, first)
first.getmm()
first.__dict__
second.getmm()
stats.get_mm = MethodType(get_mean, stats)
first.__dict__
first.get_mm()
first.getmm()
конструктор
class stats:
def __init__( self ):# Конструктор
print('first call!')
self.mean = 0
self.std = 0
def calc( self, x ):
self.mean = np.mean( x )
self.std = np.std( x )
def pri( self ):
print( self.mean, self.std)
forth = stats()
forth.pri()
class stats:
def __init__( self, x = None ):# Конструкотор, расчитан на два случая
if x is None:
print('simple call!')
self.mean = 0
self.std = 0
return
print('data call!')
self.calc( x )
def calc( self, x ):
self.mean = np.mean( x )
self.std = np.std( x )
def pri( self ):
print( self.mean, self.std)
fifth = stats()
fifth.pri()
sixth = stats( a )
sixth.pri()
Переменная класса
Важно, что ранее объявленные переменные были имено частью self, а не класса. Иначе будет следующий эффект
class stats:
all_means = []
def __init__( self, x = None ):# Конструкотор, расчитан на два случая
if x is None:
print('simple call!')
self.mean = 0
self.std = 0
return
print('data call!')
self.calc( x )
def calc( self, x ):
self.mean = np.mean( x )
self.std = np.std( x )
self.all_means.append( self.mean )
def pri( self ):
print( self.mean, self.std)
first = stats( a )
second = stats( b )
print( first.all_means, second.all_means ) # Переменная all_means одна и таже для всех объектов.
Наследование
class collect( stats ):
def __init__(self):
self.data = []
def add( self, x ):
self.data.append( x )
def calc_stats( self ):
self.calc( self.data )
#super( collect, self).calc( self.data )#Можно и так.
acum = collect()
acum.add( 1. )
acum.add( 2. )
acum.data
acum.calc_stats()
acum.pri()
acum.add( -2. )
acum.calc_stats()
acum.pri()
Подмена вызова
class stats:
all_means = []
def __init__( self, x = None ):# Конструкотор, расчитан на два случая
if x is None:
print('simple call!')
self.mean = 0
self.std = 0
return
print('data call!')
self.calc( x )
def filt( self, x):
return x
def calc( self, x ):
x = self.filt( x ) # Использует метод filt
self.mean = np.mean( x )
self.std = np.std( x )
self.all_means.append( self.mean )
def pri( self ):
print( self.mean, self.std)
class collect( stats ):
def __init__(self):
self.data = []
def filt( self, x): # Новый фильтр. Он подменяет из родительского.
xx = [ qq for qq in x if qq > 0]
return xx
def add( self, x ):
self.data.append( x )
def calc_stats( self ):
self.calc( self.data )
acum_pos = collect()
acum_pos.add( 1. )
acum_pos.add( 2. )
acum_pos.add( -2. )
acum_pos.calc_stats() # При подсчете статистик будет вызван новый фильтр.
acum_pos.pri() # Учитывались только положительные!
x = models.norm.rvs(2, 3, size = 100)
models.norm.fit( x )
np.mean( x )
np.std( x ) # Вот и все, т.е. по формулам.
pi = 0.3
lmbd = 3
def zip_pmf( x, pi = pi, lmbd = lmbd ):
if pi < 0 or pi > 1 or lmbd <= 0:
return np.zeros_like( x )
return (x == 0) * pi + (1 - pi) * models.poisson.pmf( x, lmbd)
zip_pmf( 0 )
zip_pmf( 1 )
models.poisson.pmf( 0, lmbd)
fig, ax = plt.subplots(figsize=(8, 6))
xs = np.arange(0, 10)
palette = ['g','b']
ax.bar( 2.5*xs, models.poisson.pmf( xs, lmbd), width = 0.9, color = palette[0], label='Poisson')
ax.bar( 2.5*xs + 1, zip_pmf( xs ), width = 0.9, color = palette[1], label='Zero - Poisson')
ax.set_xticks( 2.5 * xs + 1 -0.6*0.9)
ax.set_xticklabels( xs )
ax.set_xlabel('$x$')
ax.set_ylabel('$P(X=x)$')
ax.legend()
N = 1000
iz = models.bernoulli.rvs( pi, size = N)
x = (1 - iz) * models.poisson.rvs( lmbd, size = N)
print(x[::50])
x.shape
fig, ax = plt.subplots(figsize=(8, 6))
ax.hist( x, width=0.8, bins = np.arange( x.max() + 1 ), normed = True)
ax.set_xticks( np.arange( x.max() + 1 ) + 0.4 )
ax.set_xticklabels( np.arange( x.max() + 1 ) )
ax.set_xlabel('$x$')
ax.set_ylabel('Частота')
from statsmodels.base.model import GenericLikelihoodModel
class ZeroInflated(GenericLikelihoodModel):
def __init__(self, endog, exog = None, **kwds):
if exog == None:
exog = np.zeros_like( endog )
print("init !!")
super( ZeroInflated, self).__init__( endog, exog, **kwds)
def nloglikeobs( self, params):
pi = params[0]
lmbd = params[1]
value = -np.log( zip_pmf( self.endog, pi = pi, lmbd = lmbd) )
print("log, ", np.sum(value))
return value
def fit(self, start_params = None, maxiter = 10000, maxfun = 5000, **kwds):
if start_params is None:
lmbd_start = self.endog.mean()
extra_zeros = ( self.endog == 0 ).mean() - models.poisson.pmf(0, lmbd_start )
start_params = np.array([extra_zeros, lmbd_start])
return super( ZeroInflated, self ).fit(start_params = start_params, maxiter = maxiter, maxfun = maxfun, **kwds)
model = ZeroInflated( x )
#model.fit()
ZeroInflated.fit(model)
res = model.fit()
res.params # Profit!
Используем распределение Бернули
Задается вероятностью p происхождения события. Например, вероятность кликнуть на банер, зайти на сайт, прийти в магазин или купить товар. Это все разные вероятности. Важное, что некое событие может произойти и неоднократно.
models.bernoulli.rvs( 0.7, 10) # Так нельзя, так как 10 будет использоватся как сдвиг.
В питоне функции могут использовать два типа аргументов или иначе параметров, которые ей передаются. Первый тип позиционый, который является стандартным для большинства языков программиования. Суть в том, что параметры присваиваются в порядке и задания при вызове. А бывают аргументы именные, которые можно присвоить вне очереди.
Рассмотреть документацию по данной функции https://docs.scipy.org/doc/scipy/reference/generated/scipy.stats.bernoulli.html
Конкретно изучим описание функции rvs: rvs(p, loc=0, size=1, random_state=None)
Последнее к сожалению позволит понять какие значения были присвоены аргументам по умолчанию. часть аргументов просто перечислена (а данном случае только p), а часть через знак '='.
# Вообще у любого объекта, функции по соглашению может быть поле __doc__ которое дает краткую документацию по нему.
# Так как это поле является строчкой, то для аккуратного вывода лучше применить функцию print
# Подробнее о ней будет позже. Сейчас она используется для вывода самоописания функции.
print( models.bernoulli.rvs.__doc__ )
Из данного описания следует, что сначала идут позиционные параметры (arg1, arg2, ...), а потом именные (loc, size, random_state). Тем не менее, если они идут подряд, то присвоение идет подряд.
# В данном случае:
models.bernoulli.rvs( 0.7, 10) # 0.7 это параметр распределения, а 10 присвоилось аргументу loc.
models.bernoulli.rvs( 0.7, 10, 5) # Вот теперь очередь дошла до размера, которому было присвоено значение 5.
# Но можно было бы и явно присвоить переменной отвечащей за размер размер:
models.bernoulli.rvs( 0.7, size = 10 ) # Аргументу size присвоено значение вне очереди.
# Нельзя позиционный аргумент подавать после именных. Последнее выдаст соответствующую ошибку:
models.bernoulli.rvs( size = 10, 0.7 )
# Но если аргументы заданы по имени, то их можно подавать функции в произвольном порядке.
models.bernoulli.rvs( size = 10, p = 0.7 )
models.bernoulli.rvs( size = 10, p = 0.7, loc = 5 ) # В данном случае было задано и значение параметра loc, положения.
Закон больших чисел
Нам нужно взять среднее количество произошедших событий. Точнее понять как это среднее рапределено. Теория говрит, что оно должно стремится к нормальному распределению.
# Имитация закона больших чисел:
data = [] # Содали пустой список
for i in range( 1000 ): # Мы явно указали размер нешей выборки
lap = models.bernoulli.rvs( 0.7, size=100 ) # а также количество элементов в устреднение.
data.append( np.sum( lap )/100 )
plt.hist( data );
# Хотим цикл ещё раз запустить с другим количеством эллементов в сумме.
Продолжаем с Бернули
# Вводим понятие функции, т.е. сами пишем.
def buildBernSum( p, cycles = 1000, elems = 100 ): # Два аргумента. Оба необязательны. Имеют значение по умолчанию.
data = [] # Содали пустой список
for i in range( cycles ):
lap = models.bernoulli.rvs( p, size=elems ) # 0.7
data.append( np.sum( lap )/elems )
return data # Вернули результат.
bern = buildBernSum( 0.7 )
plt.hist( bern );
bern = buildBernSum( 0.7, elems = 10000 )
plt.hist( bern );
bern = buildBernSum( elems = 10000, p = 0.7)
plt.hist( bern );
bern = buildBernSum( 0.7, 100000, 10000 ) # Позиционное присвоение параметров
plt.hist( bern );
bern = buildBernSum( 0.7, 4000, 3000 )
plt.hist( bern );
#С увеличением параметров распределение будет стремится к нормальному
Anderson-Darling
Предыдущая статистика работает только для нормального распределения. Для данного же теста можно задать используемое распределение: norm, expon, logistic, gumbel, gumbel_l, gumbel_r, extreme1. Пока будем оставатся в рамках нормального.
В Шапиро-Уилка граница была единственной. Той что обычно принято пользоваться. В методе же Андерсона таких границ считается несколько для различных уровней значимости.
dat = models.uniform.rvs( size = 100 ) # Пвторно проверим на равномером распределении.
plt.hist( dat );
models.anderson( dat, dist='norm' ) # Тест отвергает данные как принадлежащие равномерному распределению.
# Статистику нужно возвести в квадрат и сравнить с критерием. Если больше, то гипотеза отвергается.
Проверка того что распределение то что нужно.
dat = models.norm.rvs( 2, 3, 100 ) # Для начала прверим для нормального распределения.
res = models.anderson( dat )
res
res[1] # Значение значимости критических уровней.
res[2] # Сами уровни.
# Можно представить их матрицей.
table = np.array([res[1], res[2]])
table
# Хотим красивое представление результата/критерия в виде таблицы
import pandas as pd # Очередная библиотека полезная для представления данных в таблицах, и их считывания.
pd.DataFrame( table, columns=['Уровни', 'Значимость'] ) # Система ругнется так как матрица не правильного размера.
# Данная функция ожидала матрицу состоящую из двух столбцов.
# В Numpy естественно есть функция для транспонирования матрицы.
np.transpose( table )
# Транспонироваие дает результат и в общем случае: она переставляет оси.
x = np.ones((1, 2, 3)) # ones вспомогательная функция составляющая многомерную матрицу из единиц. В скобках указан размер.
x
xT = np.transpose(x, (1, 2, 0)) # В скобках указана перестановка исходных осей, т.е. оси будут идти в порядке 1, 2, 0.
xT
x.shape, xT.shape # Если оси переставлены, то мы это увидем на соответствующем образом изменившимся порядке их размере.
# Вначале должна идти первая компонента, потом вторая, а уже потом нулевая векто рараеров (1, 2, 3).
# Тогда получим 2, 3, 1
# Тогда:
pd.DataFrame( np.transpose( table ), columns=['Значение статистики', 'Уровень значимости'] )
# Получаем красивую таблицу.
# Но можно иначе.
table = np.dstack( (res[1], res[2]) )
table, table.shape # Но появится лишняя размерность.
# Её можно убрать взяв первый элемент:
table[0], table[0].shape
# тогда опять можем получить таблицу.
pd.DataFrame( table[0], columns=['Значение статистики', 'Уровень значимости'] )
# Скорее всего ошибка в коде. Значимости нужны в обратном порядке.
table = np.dstack( (res[1], res[2][::-1] ) )[0]
pd.DataFrame( table, columns=['Значение статистики', 'Уровень значимости'] )
# В res[0] у нас значение статистики.
res[0] # Для проверки гипотезы нужно сравнить это число с критическими значениями.
Бернули
bern = buildBernSum( 0.7, 4000, 3000 )
res = models.anderson( bern ) # Второй параметр указывает на распределение. Поумолчанию проверка на нормальность.
res # Дествительно стремимся к нормальному распределению.
params = models.norm.fit( bern )
params
models.norm.interval( 0.95, params[0], params[1] )
models.bernoulli.interval( 0.7, 0.7 ) # Системе плохо.
models.norm.interval( 0.95 )
pp = np.mean( bern )
pp
dd = np.sqrt(pp*(1-pp)/3000)
dd
lim = models.norm.interval( 0.95, 0, 1 )
lim
pp + dd*lim[0], pp + dd*lim[1]
models.laplace.rvs()
import statsmodels.api as sm
import statsmodels.sandbox.regression.predstd as predstd
nsample = 100
x = np.linspace(-5, 5, 100)
X = np.column_stack((x, x, x**2)) # Отмечу, что переменная x повторена. т.е. есть избыток в данных.
beta = np.array([1, 0.5, 0.5, 0.5])
e = np.random.normal(size=nsample)
#X
sm.
X = sm.add_constant(X)
y = np.dot(X, beta) + e
model = sm.OLS(y, X)
results = model.fit()
print(results.summary())
print('Parameters: ', results.params)
print('R2: ', results.rsquared)
model = sm.OLS(y, X)
results = model.fit()
print(results.summary())
yy = np.dot( X, results.params)
plt.plot( x, yy)
plt.plot( x, y, ".")
results = sm.OLS(y, sm.add_constant(x)).fit() # Все внимание на F-statistic, Porb(...), R-squared.
print(results.summary()) # Они показывают, что мы плохо приблизили. "Слуйчайная прямая" дала бы схожий результат.
yy = np.dot( sm.add_constant(x), results.params)
plt.plot( x, yy)
plt.plot( x, y, ".")
Доверительный интервал
model = sm.OLS(y, X) # Возвращаеся к прошлой модели.
results = model.fit()
yy = np.dot( X, results.params)
prstd, iv_l, iv_u = predstd.wls_prediction_std( results )
plt.plot(x, yy, 'g-', label="data")
plt.plot(x, y, 'b.', label="True")
plt.plot(x, results.fittedvalues, 'r--.', label="OLS")
plt.plot(x, iv_u, 'r--')
plt.plot(x, iv_l, 'r--')
plt.legend(loc='best');
R = [[0, 1, 0, 0], [0, 0, 1, 0]]
print(np.array(R))
print(results.f_test(R)) # т.е. нестоит отбрасывать гипотезу о том, что 2 и 3 параметр идентичны.
R = [[0, 0, 0, 1], [0, 1, 0, 0]]
print(np.array(R))
print(results.f_test(R)) # первая переменная и предпоследняя явно разные.