import numpy as np # Librería de manejo de matrices y arreglos de números
from numpy.random import binomial # Generador aleatorio de números basado en la distribución binomial
from scipy.stats import binom # Permite implementar la función binomial
from math import factorial # Permite implementar la función factorial
import matplotlib.pyplot as plt # Permite realizar graficos de nuestros datos
# Definición de la deistribución binomial
def my_binomial(k, n, p):
return factorial(n)/(factorial(k) * factorial(n-k)) * pow(p, k) * pow(1-p, n-k)
print(f'My binomial: {my_binomial(2, 3, 0.5)}')
dist = binom(3, 0.5)
dist.pmf(2)
# Función de distribución acumulada para los números menores e igual a 2
print(7/8)
dist.cdf(2)
# Simulación con 100 lanzamientos de moneda equilibrada
# (ejecuta esta celda varias veces para observar la variación en los resultados)
p = 0.5
n = 3
binomial(n, p)
# Vamos a hacer un experimento generando una muestra de conjuntos de lanzamientos de a 3 monedas
arr = []
for _ in range(100):
arr.append(binomial(n, p))
def plot_hist(num_trials):
values = [0, 1, 2, 3]
arr = []
for _ in range(num_trials):
arr.append(binomial(n, p))
simulated_distribution = np.unique(arr, return_counts=True)[1]/len(arr)
teoric_distribution = [binom(3, 0.5).pmf(k) for k in values]
plt.bar(values, simulated_distribution, color='springgreen')
plt.bar(values, teoric_distribution, alpha=0.5, color='darkblue')
plt.title(f'{num_trials} experimentos')
plt.show()
print(f'Distribución Simulada: {np.unique(arr, return_counts=True)[1]/len(arr)}')
print(f'Distribución Teórica: {[binom(3, 0.5).pmf(k) for k in values]}')
plot_hist(20)
plot_hist(200)
plot_hist(20000)
import pandas as pd
from scipy.stats import norm
def gaussian(x, mu, sigma):
return 1/(sigma * np.sqrt(2 * np.pi)) * np.exp(-0.5 * pow((x-mu)/sigma, 2))
x = np.arange(-4, 4, 0.1)
y = gaussian(x, 1.0, 0.5)
plt.plot(x, y)
distribution = norm(0, 1)
x = np.arange(-4, 4, 0.1)
y = [distribution.pdf(value) for value in x] # pdf (probability density function)
plt.plot(x, y)
# Scipy nos ayuda a sacar la Distribución de probabilidad acumulada
# ya que algunas funciones no son fáciles de integrar
distribution = norm(0, 1)
x = np.arange(-4, 4, 0.1)
y = [distribution.cdf(value) for value in x] # cdf (cumulative distribution function)
plt.plot(x, y)
# Para cargar archivos excel pega el siguiente comando: pip install xlrd
# al final del archivo init.ipynb de la sección Environment
df = pd.read_excel('s057.xls')
df
# Trabajamos solo con la columna de tamaño de alas de mosca pero a partir de la posición 3
# ya que las primeras son solo etiquetas como poder en la tabla superior
arr = df['Normally Distributed Housefly Wing Lengths'].values[3:]
values, dist = np.unique(arr, return_counts=True)
plt.bar(values, dist)
# Estimación de distribución
mu = arr.mean()
sigma = arr.std()
x = np.arange(30, 60, 0.1)
dist = norm(mu, sigma)
y = [dist.pdf(value) for value in x]
plt.plot(x, y)
# Como la distribución está aplastada por ser probabilidades (van de 0 a 1)
# tenemos que normalizar nuestor conteo para mostrar las 2 gráficas juntas
values, dist = np.unique(arr, return_counts=True)
plt.bar(values, dist/len(arr))
from numpy.random import normal # Generador de números aleatorios basado en la dist. normal
sample = normal(size=10000) # generador
plt.hist(sample, bins=100)
plt.show()
sample = normal(loc=50, scale=5, size=10000) # mu = 50, sigma = 5
mu = sample.mean()
sigma = sample.std()
dist = norm(mu, sigma)
values = [value for value in range(30, 70)]
probabilities = [dist.pdf(value) for value in values]
plt.hist(sample, bins=100, density=True)
plt.plot(values, probabilities)
plt.show()
from numpy import hstack
from sklearn.neighbors import KernelDensity
# Contruimos una distribución bimodal
sample1 = normal(loc=20, scale=5, size=300)
sample2 = normal(loc=40, scale=5, size=700)
sample = hstack((sample1, sample2))
model = KernelDensity(bandwidth=2, kernel='gaussian')
sample = sample.reshape((len(sample), 1))
model.fit(sample)
values = np.asarray([value for value in range(1, 60)])
values = values.reshape((len(values), 1))
probabilities = model.score_samples(values) # probabilidad logarítmica
probabilities = np.exp(probabilities) # inversión de probabilidad
plt.hist(sample, bins=50, density=True)
plt.plot(values, probabilities)
plt.show()
from mpl_toolkits.mplot3d import Axes3D # Permite construir un gráfico 3D
from matplotlib import cm # Librería para acceder a diferentes paletas de colores
def likelihood(y, yp):
return yp*y + (1-yp)*(1-y)
fig = plt.figure()
ax = fig.gca(projection='3d')
Y = np.arange(0, 1, 0.01)
YP = np.arange(0, 1, 0.01)
Y, YP = np.meshgrid(Y, YP) # Generamos una malla de valores
Z = likelihood(Y, YP)
# Ya que la gráfica de una función de 2 variables es una superficie
surf = ax.plot_surface(Y, YP, Z, cmap=cm.coolwarm) # Creamos una superficie
fig.colorbar(surf, shrink=0.5, aspect=5) # Barra de colores para identificar los máximos
plt.show()
from sklearn.datasets import load_iris
from sklearn.linear_model import LogisticRegression
atrib_names = ['sepal length', 'sepal width', 'petal length', 'petal width']
X, y = load_iris(return_X_y=True)
X[:10]
y[:100] # Como queremos un problema binario solo escogemos las 2 primeras clases 0 y 1
clasif = LogisticRegression(random_state=10, solver='liblinear').fit(X[:100], y[:100])
# random_state: generador aleatorio que inicializa las variables de nuestro modelo
# solver: elegimos el método de optimización de nuestro modelo, en este caso *liblinear*
clasif.coef_
model_coefs = pd.DataFrame(clasif.coef_, columns=atrib_names)
model_coefs
clasif.predict(X[:100, :])
clasif.predict_proba(X[:10, :])
from math import factorial
# Definición de la distribución binomial
def my_binomial(k, n, p):
return factorial(n)/(factorial(k)*(factorial(n-k)))*pow(p,k)*pow(1-p, n-k)
def my_cumulative_binomial(k, n, p):
total = 0
for i in range(k+1):
total += my_binomial(i, n, p)
return total
my_binomial(3, 12, 0.5)
my_cumulative_binomial(5, 10, 0.5)
my_cumulative_binomial(5, 10, 0.5)
my_binomial(3, 12, 0.3)
my_cumulative_binomial(5, 10, 0.3)
my_cumulative_binomial(5, 10, 0.3)
import numpy as np
from scipy.stats import norm
import matplotlib.pyplot as plt
mu = 4
sigma = 0.1
distribution = norm(mu, sigma)
x = 4
y = distribution.pdf(x)
print(y)
x = -10
y = distribution.pdf(x)
print(y)
x = 10
y = distribution.pdf(x)
print(y)
X = np.arange(3.5, 4, 0.01)
y = [distribution.cdf(value) for value in X]
plt.plot(X, y)
print('P(X <= 4) =',distribution.cdf(4))
X = np.arange(4, 4.5, 0.01)
y = [distribution.cdf(value) for value in X]
plt.plot(X, y)
print('P(X >= 4) =', round(distribution.cdf(5)- distribution.cdf(3.9999), 2))
import numpy as np
from numpy.random import binomial
import matplotlib.pyplot as plt
def generate_binomial_trial(trials, coin_toss):
'''
El resultado de esta función es generar un conjunto de
experimentos binomiales (trials) y de cada uno obtener las
cantidades de éxitos en cada secuencia (ej. lanzar monedas).
* trial: es una secuencia de <coin_toss> lanzaientos de moneda
* coin_toss: es el número de monedas lanzadas en cada trial
'''
arr = []
for _ in range(trials):
arr.append(binomial(coin_toss, 0.5))
values, dist = np.unique(arr, return_counts=True)
return arr, values, dist
arr, values, dist = generate_binomial_trial(100000, 100)
plt.bar(values, dist)
mu = np.mean(arr)
sigma = np.std(arr)
distribution = norm(mu, sigma)
x = np.arange(20, 80, 0.1)
y = [distribution.pdf(value) for value in x]
plt.plot(x, y, color='orange')
# Datos
plt.bar(values, dist/len(arr))
import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
# Recuerda que este archivo lo puedes bajar de:
# https://seattlecentral.edu/qelp/sets/057/057.html
df = pd.read_excel('s057.xls')
arr = df['Normally Distributed Housefly Wing Lengths'].values[3:]
values, dist = np.unique(arr, return_counts=True)
print(values)
plt.bar(values, dist)
def optimal_mu(arr=arr):
mu = 0
for i in arr:
mu += i
return mu/len(arr)
def optimal_sigma(arr=arr):
mu = optimal_mu(arr)
sigma = 0
for i in arr:
sigma += pow(i - mu, 2)
return np.sqrt(sigma/len(arr))
print(optimal_mu(), arr.mean())
print(optimal_sigma(), arr.std())
from scipy.stats import norm
values, dist = np.unique(arr, return_counts=True)
plt.bar(values, dist/len(arr))
dist = norm(optimal_mu(), optimal_sigma())
x = np.arange(30, 60, 0.1)
y = [dist.pdf(value) for value in x]
plt.plot(x, y, color='orange')
from sklearn.datasets import make_blobs
import matplotlib.pyplot as plt
import numpy as np
from scipy.stats import norm
X, y = make_blobs(n_samples=10000, centers=2, n_features=2, random_state=1)
# Esta función ajusta una gaussiana
# a un conjunto de data
def fit_distribution(data):
mu = data.mean()
sigma = data.std()
dist = norm(mu, sigma)
return dist
plt.scatter(X[y==1][:, 0], X[y==1][:,1], label='1', color='red')
plt.scatter(X[y==0][:, 0], X[y==0][:,1], label='0', color='blue')
plt.legend()
# Calculamos priors
def prior(c):
return len(X[y==c])/len(X)
# Tenemos 4 posibles distribuciones a ajustar (verosimilitud)
def dist_X0(c):
if c == 0:
return fit_distribution(X[y==0][:, 0])
elif c == 1:
return fit_distribution(X[y==1][:, 0])
def dist_X1(c):
if c == 0:
return fit_distribution(X[y==0][:, 1])
elif c == 1:
return fit_distribution(X[y==1][:, 1])
def likelihood(X, c):
return dist_X0(c).pdf(X[0]) * dist_X1(c).pdf(X[1])
# Posterior
def probability(c, X):
return prior(c) * likelihood(X, c)
predictions = [np.argmax([probability(0, vector), probability(1, vector)]) for vector in X]
from sklearn.metrics import confusion_matrix
confusion_matrix(y, predictions)
def class_distribution(x, y):
print()
return np.argmax([probability(0, [x, y]), probability(1, [x, y])])
class_distribution(-6, 0)
class_distribution(-4, 0)
plt.scatter(X[y==1][:, 0], X[y==1][:, 1], label='1', color='red', marker='*')
plt.scatter(X[y==0][:, 0], X[y==0][:, 1], label='0', color='blue', marker='*')
plt.scatter(-6, 0, color='red', marker='s', s=53)
plt.scatter(-4, 0, color='blue', marker='s', s=53)
plt.legend()