# Dependencias
import numpy as np
from numpy.random import binomial
from scipy.stats import binom
import scipy.stats
from math import factorial
import matplotlib.pyplot as plt
def my_binomial(k,n,p):
return factorial(n)/(factorial(k)*factorial(n-k))*pow(p,k)*pow(1-p,n-k)
my_binomial(2,3,0.5)
# scipy.stats.binom(numero de intentos, probabilidad).pmf(numero de exitos)
dist = binom(3, 0.5)
# pmf = probability mass function = funcion de densidad de probabilidad
pmf = dist.pmf(2)
probabilidades_individuales = scipy.stats.binom(3,0.5).pmf(range(0,3))
probabilidades_sumadas = round(np.sum(probabilidades_individuales),3)
print(f"""
Probabilidades Individuales: {probabilidades_individuales}
Probabilidades Sumadas: {round(probabilidades_sumadas,3)}
""")
# scipy.stats.binom(numero de intentos = 3 , probabilidad de eventos = 0.5).cdf(casos exitoso = 2)
# Cumulative distribution function. = funcion de distribucion acumulada
dist.cdf(2)
# Comprobamos que los resultados sean correctos
print(f"""
{probabilidades_sumadas}
{dist.cdf(2)}
{7/8}
""")
# simulacion con 100 lanzamientos de moneda equilibrada
# (ejecuta esta celda varias veces para observar la variacion de los resultados)
p = 0.5
n = 3
binomial(n,p) #numpy.random.binomial
# Resultados de probabilidades teóricas, mientras más intentos hagamos, más se acerca a las probabilidades de la escuela frecuentista
arr = []
def simulada(intentos,exito,probabilidad):
for i in range(intentos):
arr.append(binomial(exito, probabilidad))
x,y = np.unique(arr, return_counts=True)
plt.bar(x,y)
plt.show()
print(np.unique(arr, return_counts=True)[1]/len(arr)) # imprimimos las probabilidades
simulada(50000,3,.5) # Si subimos el numero de intentos, las probabilidades quedaran tal como la distribucion acumulada de probabilidades individuales
values=[0,1,2,3]
[binom(3,.5).pmf(k) for k in values]
arr = []
def simulada(intentos, exito, probabilidad):
for _ in range(intentos):
arr.append(binomial(exito, probabilidad))
x = np.copy(values)
y = [binom(3,.5).pmf(k) for k in values]
plt.bar(x,y)
plt.show()
print(y) # imprimimos las probabilidades
simulada(10000,3,.5)
def plot_hist(num_trials):
values = [0,1,2,3]
arr = []
for _ in range(num_trials):
arr.append(binomial(3,0.5))
distribucion_simulada = np.unique(arr, return_counts=True)[1]/len(arr)
distribucion_teorica = [binom(3,0.5).pmf(k) for k in values]
plt.bar(values, distribucion_simulada,alpha=0.5, color = 'red')
plt.bar(values, distribucion_teorica, alpha=0.5,color='green')
plt.title('{} experimentos'.format(num_trials))
plt.show()
plot_hist(20)
plot_hist(200)
plot_hist(200000)
# Dependencias
import pandas as pd
import numpy as np
import matplotlib.pylab as plt
from scipy.stats import norm
def gaussiana(x, median, std):
return 1/(std*np.sqrt(2*np.pi))*np.exp(-0.5*pow((x-median)/std,2))
x = np.arange(-4,4,0.1)
y = gaussiana(x, 0.0, 1.0)
plt.plot(x,y)
x = np.arange(-4,4,0.1)
y1 = gaussiana(x, -2.0, 1.0) # media de -2
y2 = gaussiana(x, -1.0, 1.0) # media de -1
y3 = gaussiana(x, 1.0, 1.0) # media de 1
y4 = gaussiana(x, 2.0, 1.0) # media de 2
plt.plot(x,y1)
plt.show()
plt.plot(x,y2)
plt.show()
plt.plot(x,y3)
plt.show()
plt.plot(x,y4)
plt.show()
x = np.arange(-4,4,0.1)
y1 = gaussiana(x, 0, 2.0) # Desviacion estandar de 2
y2 = gaussiana(x, 0, 1.0) # Desviacion estandar de -2
y3 = gaussiana(x, 0, 0.5) # Desviacion estandar de -1
y4 = gaussiana(x, 0, 0.1) # Desviacion estandar de 1
plt.plot(x,y1)
plt.show()
plt.plot(x,y2)
plt.show()
plt.plot(x,y3)
plt.show()
plt.plot(x,y4)
plt.show()
# scipy.stats.norm(mean, std)
# pdf = Probability density function.
dist = norm(0,1)
x = np.arange(-4,4, 0.1)
y = [dist.pdf(value) for value in x]
plt.plot(x,y)
# cdf = Cumulative distribution function.
dist = norm(0,1)
x = np.arange(-4,4,0.1)
y = [dist.cdf(value) for value in x]
plt.plot(x,y)
# Importamos la base de datos
df = pd.read_excel('./data/s057.xls')
df.head()
# Guardamos solo la columna de la longitud de alas que tiene todos los datos, omitiendo los 3 primeros valores que no nos sirven para la distribución normal
arr = df['Normally Distributed Housefly Wing Lengths'].values[3:]
print(arr)
# Separamos los valores únicos del array con NumPy y activamos que cuente cuantas veces se repite cada valor único
values, dist = np.unique(arr, return_counts=True)
print(f"""
{values} Lista de valores único
{dist} Lista de frecuencia de valores unicos
""")
# Graficamos con matplotlib
plt.bar(values,dist)
plt.xlabel('Largo de alas en milímetros')
plt.ylabel('frecuencia')
# El gráfico que obtuvimos a simple vista luce como una campana de gauss, pero para asegurarnos de que esto sea así, lo verificaremos haciendo una estimación de una distribución
mean = arr.mean() # Calculamos el promedio del array de datos
std = arr.std() # Calculamos la desviación estándar de los datos
x = np.arange(30, 60, 0.1) # Creamos una lista de valores que se encuentre un poco mas alejado de los valores limites de la distribución (30)<-35 || 55->(60)
dist = norm(mean, std) # Con scipy definimos la distribución normal con el promedio y desviacion estandar de los datos
y = [dist.pdf(value) for value in x] # Calculamos Y con la densidad de probabilidad en cuanto a los datos de dist
plt.plot(x,y) # Graficamos la estimación de la distribución
values, dist = np.unique(arr, return_counts=True)
plt.bar(values,dist/len(arr)) # Graficamos la distribución del ejercicio anterior, y normalizamos la lista dist con el total de espacios en el array original de los datos
plt.show()
# Dependencias
import numpy as np
from numpy.random import normal
from scipy.stats import norm
import matplotlib.pyplot as plt
sample = normal(size = 10000) # generador aleatorio
plt.hist(sample, bins = 30);
sample = normal(loc=50, scale=5, size=10000) # mean=50, std=5, size=10000
plt.hist(sample,bins=30, density=True);
mean = sample.mean()
std = sample.std()
dist = norm(mean,std)
values = [value for value in range(30,70)]
probabilidades = [dist.pdf(value) for value in values]
plt.plot(values, probabilidades)
plt.show()
plt.hist(sample, bins=30, density=True)
plt.plot(values, probabilidades)
plt.show()
# Codigo completo
sample = normal(loc=50, scale=5, size=5000000) # 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=50, density=True)
plt.plot(values, probabilities)
plt.show();
# Dependencias
import matplotlib.pyplot as plt
from numpy import hstack
from sklearn.neighbors import KernelDensity
# construccion de una distribucion bimodal
sample1 = normal(loc=20, scale=5, size=300)
sample2 = normal(loc=40, scale=5, size=700)
sample = hstack((sample1, sample2))
# graficamos
plt.hist(sample, bins=100);
model = KernelDensity(bandwidth=2, kernel='gaussian') # (parametro de suavizado, funcion base)
# print(sample)
sample = sample.reshape((len(sample), 1))
# print(sample)
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 para optimizar calculos computacionales
probabilities = np.exp(probabilities) # invertimos los resultados de las probabilidades logaritmicas con su exponencial para tener la escala original
plt.hist(sample, bins=50, density=True)
plt.plot(values, probabilities)
plt.show()
# Dependencias
from scipy.stats import norm
import numpy as np
import matplotlib.pyplot as plt
from sklearn.neighbors import KernelDensity
from datetime import date, timedelta
import requests
import json
# Valores de entrada
start_date = date.today() - timedelta(days=7)
end_date = date.today()
api_key = 'JKsQ8eXx4qD051uWOX4Q23ptCbGzr0B7e5DpcDKp'
api_url = 'https://api.nasa.gov/neo/rest/v1/feed?start_date={start_date}&end_date={end_date}&api_key={api_key}'.format(start_date = start_date, end_date = end_date, api_key = api_key)
# Consumo de API
def obtain_asteroid_speed():
try:
request = requests.get(api_url)
data = json.loads(request.text)
except:
print('xD')
speeds = []
for day in data['near_earth_objects']:
for obj in data['near_earth_objects'][day]:
speeds.append(round(float(obj['close_approach_data'][0]['relative_velocity']['kilometers_per_second'])))
speeds = np.array(speeds)
return speeds
arr = obtain_asteroid_speed()
# Forma de los datos
plt.hist(arr, bins=25);
def parametric_estimation_asteroid_speed():
mean = arr.mean()
std = arr.std()
x = np.arange(0, 50, 0.1)
dist = norm(mean, std)
y = [dist.pdf(value) for value in x]
values, dist = np.unique(arr, return_counts=True)
plt.bar(values,dist/len(arr))
plt.plot(x,y)
plt.xlabel('kilometer per second')
plt.ylabel('frecuency')
plt.show()
parametric_estimation_asteroid_speed()
arr = arr.reshape((len(arr),1))
def kernel_density_estimation():
model = KernelDensity(bandwidth=2, kernel='epanechnikov')
model.fit(arr)
values = np.asarray([value for value in range(1,35)])
values = values.reshape((len(values),1))
probabilities = model.score_samples(values)
probabilities = np.exp(probabilities)
plt.hist(arr, bins=25,density=True)
plt.plot(values,probabilities)
plt.xlabel('kilometer per second')
plt.ylabel('frecuency')
plt.show()
kernel_density_estimation()
# Dependencias
from mpl_toolkits.mplot3d import Axes3D
import matplotlib.pyplot as plt
from matplotlib import cm
import numpy as np
import pandas as pd
# Funcion de verosimilitud para reg. logistica
def likelihood(y,yp):
return yp*y + (1-yp)*(1-y)
fig = plt.figure()
ax = plt.axes(projection='3d')
# generamos X(vector), Y(vector) y Z(matriz)
Y = np.arange(0,1, 0.01)
YP = np.arange(0,1, 0.01)
Y, YP = np.meshgrid(Y, YP)
Z = likelihood(Y, YP)
# Superficie
surf = ax.plot_surface(Y,YP,Z, cmap=cm.coolwarm)
fig.colorbar(surf, shrink=0.5, aspect=5)
plt.show()
# La libreria que nos ayudaremos es plotly
#pip install plotly
import plotly.graph_objects as go
# generamos X(vector), Y(vector) y Z(matriz)
y = np.arange(0, 1, 0.01)
yp = np.arange(0, 1, 0.01)
Y, YP = np.meshgrid(y, yp)
Z = likelihood(Y, YP)
# graficar
fig = go.Figure(data=[go.Surface(z=Z, x=y, y=yp)])
fig.update_traces(contours_z=dict(show=True, usecolormap=True,
highlightcolor="limegreen", project_z=True))
fig.update_layout(title=' ', autosize=False,
width=500, height=500,
margin=dict(l=65, r=50, b=65, t=90))
fig.show()
# Dependencias
from sklearn.datasets import load_iris
from sklearn.linear_model import LogisticRegression
# atributos del dataset iris
atrib_names = ['sepal_length', 'sepal_width', 'petal_length', 'petal_width']
# X: atributos de los objetos
# y: tipos de flores
X, y = load_iris(return_X_y=True)
# Miramos los primeros 2 atributos del datapoint
print(X[:4])
print(y[:4])
# Instanciamos el modelo de clasificacion logistico para las dos primeras clases x: [0,1] y para los 100 primeros atributos de cada datapoint
clf = LogisticRegression(random_state=10, solver='liblinear').fit(X[:100],y[:100])
# Vemos la matriz de coeficientes
clf.coef_
# Logistic Regresion
# Preparamos los datos de entrada y salida
#X: se mantiene como esta
y = y.reshape((len(y), 1))
X[:2]
y[:2]
# Logistic Regressions
from sklearn.linear_model import LogisticRegression
# Ajustamos los datos al modelo de regresion logistica
model = LogisticRegression().fit(X,y)
# Exatitud del modelo
# Devuelve la precisión media en los datos de prueba y las etiquetas dadas.
model.score(X,y)
# Realizamos la prediccion
expected = y
predicted = model.predict(X) # lo ideal seria tener nuevos inputs, pero para este caso de estudio usaremos el mismo
predicted # Resultado de las predicciones para las flores setosa, versicolor y virginica
from sklearn import metrics
print(metrics.classification_report(expected, predicted))
import seaborn as sns
iris = load_iris()
iris = pd.DataFrame(data = np.c_[iris['data'], iris['target']], columns=atrib_names + ['species'])
iris['species'] = iris['species'].map({ 0:'setosa', 1:'versicolor', 2:'virginica'})
sns.FacetGrid(data=iris, hue='species', height=4).map(plt.scatter, 'petal_length', 'sepal_width').add_legend()
print(metrics.confusion_matrix(expected, predicted))