## Imports


```python
#Imports
import numpy as np
import matplotlib as mpl
import matplotlib.pyplot as plt
import math
import pandas as pd
import seaborn as sns
import scipy as sp
import sympy as sm
import plotly.graph_objects as go
import plotly.express as px
# from google.colab import files
from mpl_toolkits import mplot3d

from sklearn.datasets import make_blobs
from sklearn.cluster import KMeans
from sklearn.cluster import DBSCAN
from scipy.cluster.hierarchy import dendrogram, linkage, fcluster
from sklearn.preprocessing import StandardScaler

#Warnings
import warnings
warnings.filterwarnings("ignore", category=np.VisibleDeprecationWarning)

#Sinonimos
sym = sm
ip = sp.interpolate
la = np.linalg
cos = np.cos
sin = np.sin
e = np.e
pi = np.pi
cm = mpl.cm
Voronoi = sp.spatial.Voronoi
voronoi_plot_2d = sp.spatial.voronoi_plot_2d

#Styles
plt.style.use('ggplot')
plt.rcParams['figure.figsize'] = [16, 7]
plt.rcParams['axes.grid'] = True
mpl.rcParams['axes.facecolor'] = 'F2D7D5'
mpl.rcParams['axes.prop_cycle'] = mpl.cycler(color=['c', 'k', 'r'])
np.set_printoptions(precision=5, suppress=True)
```

## Guía 1
Interpolación (funciones de Tchebyshev)


```python
def T(n, a = -1, b = 1):
    '''
    Devuelve el polinomio de Chebyshev de n raíces.
    Se le puede extraer los ceros con T(n).r
    '''
    if n==0:
        return (np.poly1d([1]))
    if n==1:
        return (np.poly1d([1,0]))
    else:
        return (2*np.poly1d([1,0])* T(n-1) - T(n-2))


def cheb_zeros(n):
    '''
    Devuelve los ceros del Chebyshev de grado n como array de np
    '''
    zeros = np.polynomial.chebyshev.chebpts1(n)
    return np.array(zeros)

def cheb_zeros_plot(n):
    '''
    Devuelve los ceros de Chebyshev y los grafica con su Polinomio
    '''
    x = np.linspace(-1, 1, 200)
    zeros = np.polynomial.chebyshev.chebpts1(n)
    y = np.zeros(n)
    T = np.cos(n * np.arccos(x))
    plt.plot(zeros, y, lw=0, marker='o', ms=10)
    plt.plot(zeros, y, lw=3, color='grey')
    plt.plot(x, T, color='b')
    plt.grid()
    return np.array(zeros)

def domain_change(zeros, x_min, x_max):
    ''' Devuelve el array de zeros trasladado al nuevo dominio. zeros debe ser un np.array.
    Se puede usar directamente con np.polynomial.chebyshev.chebpts1(n) '''
    return x_min + (zeros + 1) * (x_max - x_min) / 2
```

## Guía 2
Resolución numérica de ecuaciones


```python
def Biseccion(f, a, b, n, tol=False, preview=False, p=False):
    '''
    Aproxima el valor de una raíz de f en el intervalo [a, b] utilizando el método de bisección, n veces.
    Si el ultimo intervalo es menor a tol, se devuelve el punto medio de dicho intervalo sin llegar a la iteración n
    Tanto p o preview = True van mostrando los valores para cada ciclo. Si es falso, la función devuelve sin mostrar nada.
    '''
    a = a
    b = b
    tol = tol
    for k in range(n):
        f_a = f(a)
        f_b = f(b)
        c = (a+b)/2
        f_c = f(c)
        prod = f_a * f_c
        if f_c == 0:
            return c
        if tol:
            max_err = abs(b - a)
            if max_err < tol:
                return c
        if prod > 0:
            if preview or p is True:
                print(f'Paso {k+1}: Intervalo = [{a} ; {b}]. Raíz ≈ {c}, Prod = +')
            a = c
        elif prod < 0:
            if preview or p is True:
                print(f'Paso {k+1}: Intervalo = [{a} ; {b}]. Raíz ≈ {c}, Prod = -')
            b = c
    return c

def RegulaFalsi(f, a, b, n, tol=False, preview=False, p=False):
    '''
    Aproxima el valor de una raíz de f en el intervalo [a, b] utilizando el método de bisección, n veces.
    Si el ultimo intervalo es menor a tol, se devuelve el punto medio de dicho intervalo sin llegar a la iteración n
    Tanto p o preview = True van mostrando los valores para cada ciclo. Si es falso, la función devuelve sin mostrar nada.
    '''
    a = a
    b = b
    tol = tol
    q = 0
    for k in range(n):
        f_a = f(a)
        f_b = f(b)
        c = (a * f_b - b * f_a)/(f_b-f_a)
        f_c = f(c)
        prod = f_a * f_c
        if f_c == 0:
            return c
        if tol:
            max_err = abs(b - a)
            if max_err < tol:
                return c
        if prod > 0:
            if preview or p is True:
                print(f'Paso {k+1}: Intervalo = [{a} ; {b}]. Raíz ≈ {c}, Prod = +')
            a = c
        elif prod < 0:
            if preview or p is True:
                print(f'Paso {k+1}: Intervalo = [{a} ; {b}]. Raíz ≈ {c}, Prod = -')
            b = c
        q = c
    if tol:
        print('No se alcanzó la tolerancia para ', n, 'iteraciones')
    return q

def NewtonRaphson(f, df, x0, n, tol, p=False, all_roots = False):
    tol = tol
    xn = x0
    roots = {}
    if p:
        print('x0 =', xn)
    for k in range(n):
        f_n = f(xn)
        if abs(f_n) < tol:
            if p:
                print('Raíz:', xn, '. Iteraciones:', k)
            if all_roots:
                return roots
            return xn
        df_n = df(xn)
        if df_n == 0:
            raise ValueError('df(xn) debe ser != 0')
        if all_roots:
            roots[k] = xn
        xn = xn - f_n/df_n
        if p:
            print('x'+str(k+1) + ' =', xn)
    print('No se supero la tolerancia después de ', n, 'iteraciones.')
    if all_roots:
        return roots
    return xn

def Secante(f, x0, x1, n, tol, p=False, all_roots=False):
    tol = tol
    x0 = x0
    x1 = x1
    roots = {}
    for k in range(n):
        f_0 = f(x0)
        f_1 = f(x1)
        x2 = x1 - (f_1/(f_1-f_0)*(x1-x0))
        f_2 = f(x2)
        if all_roots:
            roots[k] = x2
        if abs(f_2) < tol:
            if p:
                print('Raíz:', x2, '. Iteraciones:', k)
            if all_roots:
                return roots
            return x2
        if p:
            print(
            f'x{k}: {x0}\nx{k+1}: {x1}\nx{k+2}: {x2}\n ---------------------------------------------')
        x0 = x1
        x1 = x2

    print('No se supero la tolerancia después de ', n, 'iteraciones.')
    if all_roots:
        return roots
    return x2

def PuntoFijo(g, x0, n, tol, p=False, all_roots=False):
    x = x0
    if p:
        print('x'+str(0)+' =', x)
    if all_roots:
        all_roots = []
        all_roots.append(x)
    for k in range(n):
        x = g(x)
        if p:
            print('x'+str(k+1)+' =', x)
        if all_roots:
            all_roots.append(x)
        if abs(x-g(x)) < tol:
            if all_roots:
                return all_roots
            return x
    print('No se alcanzó la tolerancia. Último x:', x)
    if all_roots:
        return all_roots
    return x

def GaussSeidel(A, b, X0, n, tol, conv = False, p = False, all_roots = False):
    D = np.diag(np.diag(A))
    L = np.tril(A, -1)
    U = np.triu(A, 1)
    DL = D + L
    DLinv = la.inv(DL)
    X = X0
    if all_roots:
        all_roots = {}
        all_roots[0] = X
    Bgs = - DLinv @ U
    Cgs = DLinv @ b
    if conv:
        lambdas = la.eigvals(Bgs)
        if all(np.abs(l) < 1 for l in lambdas):
            print('El método de Gauss-Seidel converge para este sistema. Autovalores:\n', lambdas)
        else:
            print('Gauss-Seidel no converge para este sistema. Autovalors:\n', lambdas)
        print('--------------------------------------------------------------')
    if p:
        print('X' + str(0) + ' =\n', X)
    for k in range(1, n + 1):
        old_X = X
        X = Bgs @ old_X + Cgs
        if p:
            print('X' + str(k) + ' =\n', X)
        if all_roots:
            all_roots[k] = X
        if la.norm(X - old_X) < tol:
            if all_roots:
                return all_roots
            return X
    print('No se alcanzó la tolerancia. Último X =\n', X)
    if all_roots:
        return all_roots
    return X

def Jacobi(A, b, X0, n, tol, conv = False, p = False, all_roots = False):
    D = np.diag(np.diag(A))
    L = np.tril(A, -1)
    U = np.triu(A, 1)
    D_inv = la.inv(D)
    X = X0
    if all_roots:
        all_roots = {}
        all_roots[0] = X
    Bj = - D_inv @ (L + U)
    Cj = D_inv @ b
    if conv:
        lambdas = la.eigvals(Bj)
        if all(np.abs(l) < 1 for l in lambdas):
            print('El método de Jacobi converge para este sistema. Autovalores:\n', lambdas)
        else:
            print('Jacobi no converge para este sistema. Autovalors:\n', lambdas)
        print('--------------------------------------------------------------')
    if p:
        print('X' + str(0) + ' =\n', X)
    for k in range(1, n + 1):
        old_X = X
        X = Bj @ old_X + Cj
        if p:
            print('X' + str(k) + ' =\n', X)
        if all_roots:
            all_roots[k] = X
        if la.norm(X - old_X) < tol:
            if all_roots:
                return all_roots
            return X
    print('No se alcanzó la tolerancia. Último X =\n', X)
    if all_roots:
        return all_roots
    return X

def NR2D(F, D, x0, y0, N, tol):
    X0 = np.array([x0, y0])
    X1 = X0 - np.linalg.solve(J(X0[0], X0[1]), F(X0[0], X0[1]))
    i = 1
    print(f'i = {i}. X{i} = {X1}')
    while (np.linalg.norm(X0 - X1) > tol and i < N):
        X0 = X1
        X1 = X0 - np.linalg.solve(J(X0[0], X0[1]), F(X0[0], X0[1]))
        i = i + 1
        print(f'i = {i}. X{i} = {X1}')
    return X1
```

## Guía 3
Resolución numérica de ecuaciones difeenciales


```python
def Euler(f, x0, t0, TF, h = False, N = False):
    '''
    Sirve tanto para N como para h, pero hay que explicitarlo
    '''
    if h:
        N = int((TF - t0) / h)
    elif N:
        h = (TF - t0) / N
    else:
        raise ValueError('Debe haber valor de h o de N asignado')
    t = np.linspace(t0, TF, N+1)
    x = np.zeros(N + 1)
    x[0] = x0
    for k in range(1, N + 1):
        x[k] = x[k-1] + h * f(t[k-1], x[k-1])
    return t, x


def EulerImplicito(f, x0, t0, TF, h = False, N = False):
    '''
    Sirve tanto para N como para h, pero hay que explicitarlo
    '''
    if h:
        N = int((TF - t0) / h)
    elif N:
        h = (TF - t0) / N
    else:
        raise ValueError('Debe haber valor de h o de N asignado')
    t = np.linspace(t0, TF, N+1)
    x = np.zeros(N + 1)
    x[0] = x0
    for k in range(1, N + 1):
        x[k] = x[k-1] + h * f(t[k], x[k])
    return t, x

def Taylor2(f, ft, fx, t0, TF, x0, h = False, N = False):
    '''
    Sirve tanto para N como para h, pero hay que explicitarlo
    '''
    if h:
        N = int((TF - t0) / h)
    elif N:
        h = (TF - t0) / N
    else:
        raise ValueError('Debe haber valor de h o de N asignado')
    t = np.linspace(t0, TF, N+1)
    x = np.zeros(N + 1)
    x[0] = x0
    for i in range(1, N + 1):
        x[i] = x[i-1]+ h * f(t[i-1], x[i-1]) + (h**2 / 2) *(ft(t[i-1], x[i-1]) + fx(t[i-1], x[i-1]) * f(t[i-1], x[i-1]))
    return t, x


def RK2Especifico(f, t0, TF, x0, h = False, N = False, A1 = 0, A2 = 1, alpha = .5):
    '''
    Por defecto es Euler Modificado (A1 = 0, A2 = 1, alpha = 1)
    Sirve tanto para N como para h, pero hay que explicitarlo
    '''
    if A1 + A2 != 1:
        return ValueError('A1 + A2 debe dar exactamente 1')
    if alpha * A2 != .5:
        return ValueError('A2 * alpha debe dar exactamente 1/2')
    if h:
        N = int((TF - t0) / h)
    elif N:
        h = (TF - t0) / N
    else:
        raise ValueError('Debe haber valor de h o de N asignado')
    t = np.linspace(t0, TF, N+1)
    x = np.zeros(N + 1)
    x[0] = x0
    for i in range(1, N+1):
        x[i] = x[i-1] + h*A1*f(t[i-1], x[i-1]) + h*A2*f(t[i-1] + alpha*h, x[i-1] + alpha*h*f(t[i-1], x[i-1]))
    return t, x


def Heun(f, t0, TF, x0, h = False, N = False):
    '''
    Sirve tanto para N como para h, pero hay que explicitarlo
    '''
    if h:
        N = int((TF - t0) / h)
    elif N:
        h = (TF - t0) / N
    else:
        raise ValueError('Debe haber valor de h o de N asignado')
    t = np.linspace(t0, TF, N+1)
    x = np.zeros(N + 1)
    x[0] = x0
    for i in range(1, N+1):
        K1 = f(t[i-1], x[i-1])
        K2 = f(t[i-1] + h, x[i-1] + h*K1)
        x[i] = x[i-1] + h*(K1+K2) / 2
    return t, x

RK2 = Heun

def EulerModificado(f, t0, TF, x0, h = False, N = False):
    '''
    Sirve tanto para N como para h, pero hay que explicitarlo
    '''
    if h:
        N = int((TF - t0) / h)
    elif N:
        h = (TF - t0) / N
    else:
        raise ValueError('Debe haber valor de h o de N asignado')
    t = np.linspace(t0, TF, N+1)
    x = np.zeros(N + 1)
    x[0] = x0
    for i in range(1, N+1):
        x[i] = x[i-1] + h * f(t[i-1] + h/2, x[i-1] + h/2 * f(t[i-1], x[i-1]))
    return t, x



def RK4(f, t0, TF, x0, h = False, N = False):
    '''
    Sirve tanto para N como para h, pero hay que explicitarlo.
    Revisar las condiciones de los K
    '''
    if h:
        N = int((TF - t0) / h)
    elif N:
        h = (TF - t0) / N
    else:
        raise ValueError('Debe haber valor de h o de N asignado')
    t = np.linspace(t0, TF, N+1)
    x = np.zeros(N + 1)
    x[0] = x0
    for i in range(1, N+1):
        K1 = f(t[i-1], x[i-1])
        K2 = f(t[i-1] + h/2, x[i-1] + (h/2)*K1)
        K3 = f(t[i-1] + h/2, x[i-1] + (h/2)*K2)
        K4 = f(t[i-1] + h, x[i-1] + h * K3)
        x[i] = x[i-1] + (h/6)*(K1+2*K2 + 2*K3 + K4)
    return t, x


def EulerSistemas(F, X0, t0, TF, h = False, N = False):
    '''
    Sirve tanto para N como para h, pero hay que explicitarlo
    '''
    if len(X0) != F.__code__.co_argcount:
        raise ValueError('X0 debe tener igual cantidad de variables que F')
    if h:
        N = int((TF - t0) / h)
    elif N:
        h = (TF - t0) / N
    else:
        raise ValueError('Debe haber valor de h o de N asignado')
    n = len(X0)
    t = np.linspace(t0, TF, N+1)
    X = np.zeros((n, N + 1))
    X[:, 0] = X0
    for k in range(1, N + 1):
        X[:, k] = X[:, k-1] + h * F(t[k-1], X[:, k-1])
    return t, X

def EulerImplicitoSistemas(F, X0, t0, TF, h = False, N = False):
    '''
    Sirve tanto para N como para h, pero hay que explicitarlo
    '''
    if len(X0) != F.__code__.co_argcount:
        raise ValueError('X0 debe tener igual cantidad de variables que F')
    if h:
        N = int((TF - t0) / h)
    elif N:
        h = (TF - t0) / N
    else:
        raise ValueError('Debe haber valor de h o de N asignado')
    n = len(X0)
    t = np.linspace(t0, TF, N+1)
    X = np.zeros((n, N + 1))
    X[:, 0] = X0
    for k in range(1, N + 1):
        X[:, k] = X[:, k-1] + h * f(t[k], X[:, k])
    return t, X

def EulerModificadoSistemas(F, X0, t0, TF, h = False, N = False):
    '''
    Sirve tanto para N como para h, pero hay que explicitarlo
    '''
    if len(X0) != F.__code__.co_argcount:
        raise ValueError('X0 debe tener igual cantidad de variables que F')
    if h:
        N = int((TF - t0) / h)
    elif N:
        h = (TF - t0) / N
    else:
        raise ValueError('Debe haber valor de h o de N asignado')
    n = len(X0)
    t = np.linspace(t0, TF, N+1)
    X = np.zeros((n, N + 1))
    X[:, 0] = X0
    for i in range(1, N+1):
        X[:, i] = X[:, i-1] + h * F(t[i-1] + h/2, X[:,i-1] + h/2 * F(t[i-1], X[:,i-1]))
    return t, X

def EncontrarEquilibrios(F1, F2, xy_max = 100):
    '''
    BruteForce para encontrar equilibrios en sistemas de dos ecuaciones difs.
    Las funciones no deben estar definidas en f(x, y), NO en f(t, X)
    '''
    eqs = []
    for x in range(xy_max):
        for y in range(xy_max):
            if ((F1(x, y) == 0) and (F2(x, y) == 0)):
                eqs.append((x, y))
    return eqs
```

## Guía 4
Cadenas de Markov


```python
def MarkovIters(M, v0, n, p = True):
    nd = len(M)
    assert np.sum(M) == np.sum(np.ones(nd)), 'La matriz M no es de Markov'
    assert np.sum(v0) == np.sum(np.ones(1)), 'El vector inicial no es un vector de proporción'
    assert M.shape[0] == M.shape[1], 'La matriz M no es cuadrada'
    assert (v0.ndim == 1), 'El vector inicial no es de nx1 dimensiones (escribirlo como única fila)'
    vn = v0
    vl = 0
    if p:
        print(f'v0 = {v0}')
    for i in range(n):
        try:
            vl = vm
        except:
            pass
        vm = np.round(vn, 3)
        vn = np.round(M @ vn, 3)
        if p:
            print(f'v{i+1} = {vn}')
        if np.array_equal(vm, vn) and i >= 1:
            if p:
                print(f'v{i} y v{i+1} son iguales:\n{vm} es Estado Límite\n')
                return vm
        if np.array_equal(vn, vl):
            if p:
                print(f'v{i} y v{i+2} son iguales:\nLa Matriz no tiene Estado Límite. Alterna entre\n{vn} y {vm}')
                return None

    return vn

def MarkovEqs(M, p = True):
    nd = len(M)
    assert np.sum(M) == np.sum(np.ones(nd)), 'La matriz M no es de Markov'
    assert M.shape[0] == M.shape[1], 'La matriz M no es cuadrada'
    e = np.linalg.eig(M)
    autovals = np.round(e[0], 5)
    autovecs = e[1]
    autovals_list = list(autovals)
    autovecs = np.round((autovecs)/sum(np.abs(autovecs)), 5)
    if p:
        print('Autovalores:', autovals_list)
        print('\nAutovectores =\n', autovecs)
    eqs = []
    q = 0
    for val in autovals:
        if np.linalg.norm(val)==1:
            idx = autovals_list.index(val)
            eq = autovecs[:, idx + q]
            if p:
                print(f'\nav{q+1} = autovector para lambda {val} = \n {eq}\n')
            eqs.append(eq)
            q += 1
            autovals_list.pop(idx)
    if len(eqs) > 1:
        return eqs
    else:
        return eqs[0]
```

## Guía 5
Cuadrados Mínimos


```python
# Funcs.
def CuadradosMinimos(X, Y, R2=True, p=True):
    """Devuelve alfa y beta que minimian S(alfa, beta).
    Si R2=True, devuelve un array con [{poly1d de regresión}, {R2}]
    """
    n = len(X)
    unos = np.ones(n)
    At = np.array([unos, X])
    A = np.transpose(At)
    AtA = At @ A
    AtY = At @ Y
    params = np.linalg.solve(AtA, AtY)
    reg = np.poly1d(np.flipud(params))
    if R2:
        Y_media = sum(Y)/len(Y)
        Y_star = reg(X)
        R2 = np.linalg.norm(Y_star - Y_media, 2)**2 / np.linalg.norm(Y - Y_media, 2) ** 2
        if p:
            print(reg)
            print('R2 = ', R2)
        return np.array([reg, R2])
    if p and R2==False:
        print(reg)
        return reg
    else:
        return reg

RegresionLineal = CuadradosMinimos

def CuadradosMinimosW(X, Y, W, p=True):
    """Si R2=True, devuelve un array con [{poly1d de regresión}, {R2}]
    W debe estar en forma de vector (1D)
    """
    n = len(X)
    unos = np.ones(n)
    W = np.diag(W)
    At = np.array([unos, X])
    A = np.transpose(At)
    AtWA = At @ W @ A
    AtWY = At @ W @ Y
    params = np.linalg.solve(AtWA, AtWY)
    reg = np.poly1d(np.flipud(params))
    if p:
        print(reg)
    return reg

RegresionLinealW = CuadradosMinimosW

def CombinacionLineal(X, Y, Fs, p=True):
    lx = len(X)
    ly = len(Y)
    if lx != ly:
        raise ValueError('La cantidad de datos en x e y es distinta')
    A = np.c_[Fs[0](X)]
    for f in range(1, len(Fs)):
        A = np.c_[A, Fs[f](X)]
    At = A.T
    AtA = At @ A
    AtY = At @ Y
    alphas = np.linalg.solve(AtA, AtY)
    if p:
        for q in range(len(alphas)):
            print(f'Alfa {q} = {alphas[q]}')
    return alphas

def RegresionLinealMultiple(X1, X2, Y, R2 = False, p=True):
    lx = len(X1)
    if lx != len(X2) != len(Y):
        raise ValueError('Las dimensiones no tienen igual cantidad de datos')
    ones = np.ones(lx)
    At = np.array([ones, X1, X2])
    A = At.T
    AtA = At @ A
    AtY = At @ Y
    reg = np.linalg.solve(AtA, AtY)
    if p:
        print('Alpha 0 =', reg[0], '\nAlpha 1 =', reg[1], '\nAlpha 2 =', reg[2])
    return reg


def ObtenerR2(X, Y, regresion):
    Y_media = sum(Y)/len(Y)
    Y_star = regresion(X)
    R2 = np.linalg.norm(Y_star - Y_media, 2)**2 / np.linalg.norm(Y - Y_media, 2) ** 2
    return R2

def ObtenerR2ML(X1, X2, Y, regresion):
    Y_media = sum(Y)/len(Y)
    Y_star = regresion(X1, X2)
    R2 = np.linalg.norm(Y_star - Y_media, 2)**2 / np.linalg.norm(Y - Y_media, 2) ** 2
    return R2


def ObtenerError(X, Y, regresion):
    Y_vals = regresion(X)
    err = np.linalg.norm(Y - Y_vals, 2)
    return err**2

def ObtenerErrorML(X1, X2, Y, regresion):
    Y_vals = regresion(X1, X2)
    err = np.linalg.norm(Y - Y_vals, 2)
    return err**2
```

## Guía 6
PCA y Clustering


```python
def PCA(X_array, Return=None, theta=1.0, z_max=None, normalizar_X=True, p=False):
    '''
    Returns: Z: Todos los z\nlambdas: autovals y U\nX_star: X*\ncov o A: A\nNada: Coefs. para var objetivo
    '''
    #Calculo de X*
    l = len(X_array[:, 0])
    if normalizar_X:
        X_star =np.c_[ X_array[:, 0] - sum(X_array[:, 0])/l]
        for n in range(1, len(X_array[0])):
            col = X_array[:, n]
            X_star = np.c_[X_star, col - sum(col)/l]
        if p:
            print('X* =\n', X_star, '\n-------------------------')
    else:
        X_star = X_array

    #Calculo de Matriz de Cov
    X_star_T = X_star.T
    A = X_star_T @ X_star / l

    # Orden de autovals y autovecs
    autovals, U = np.linalg.eig(A)
    idx = autovals.argsort()[::-1]
    autovals = autovals[idx]
    U = U[:, idx]
    if p:
        print(f'Autovalores = {autovals}\nMatriz de autovectores:\n{U}\n-------------------------')
    Z = X_star @ U
    varTotal = sum(autovals)
    varAcum = autovals[0]
    varObj = theta * varTotal
    i = 1
    proporcion = 100 * varAcum / varTotal
    if p:
        print(f'Varianza total = {varTotal}')
        print(f'Varianza objetivo = {varObj}')
        print(f'Con {i} componente, la varianza es de {varAcum}\nRepresenta un {proporcion}% de la varianza total')
    while varAcum < varObj:
        varAcum += autovals[i]
        proporcion = 100 * varAcum / varTotal
        i = i+1
        if p:
            print(f'Con {i} componentes, la varianza acumulada es de {varAcum}\nRepresenta un {proporcion}% de la varianza total)')
    if p:
        print('-------------------------\nComponentes principales:\n', Z[:, 0: i])
    if Return == 'Z':
        return Z
    elif Return == 'lambdas':
        return [autovals, U]
    elif Return == 'cov' or Return =='A':
        return A
    elif Return == 'X_star':
        return X_star
    elif z_max:
        return Z[:, 0:z_max]
    else:
        return Z[:, 0:i]



def ObtenerZyXstar(X_array, p=True):
    #Calculo de X*
    l = len(X_array[:, 0])
    X_star =np.c_[ X_array[:, 0] - sum(X_array[:, 0])/l]
    for n in range(1, len(X_array[0])):
        col = X_array[:, n]
        X_star = np.c_[X_star, col - sum(col)/l]
    if p:
        print('X* =\n', X_star, '\n-------------------------')

    #Calculo de Matriz de Cov
    X_star_T = X_star.T
    A = X_star_T @ X_star / l

    # Orden de autovals y autovecs
    autovals, U = np.linalg.eig(A)
    idx = autovals.argsort()[::-1]
    autovals = autovals[idx]
    U = U[:, idx]
    Z = X_star @ U
    if p:
        print(f'Z = \n', Z)
    return np.array([Z, X_star])

def PlotZvsXstar(Z, X_star, z_min = -5, z_max = 5):
    n = 0
    m = 0
    k = 0
    x_plot = np.linspace(z_min, z_max, 100)
    for z in Z.T:
        n += 1
        m = 0
        for x in X_star.T:
            k += 1
            plt.subplot(len(X_star.T), len(Z.T), k)
            m += 1
            plt.plot(z, x, 'k.')
            plt.xlabel(f'z{n}')
            plt.ylabel(f'x{m}*', rotation=0)
            f, R2 = RegresionLineal(z, x, p=False)
            plt.plot(x_plot, f(x_plot), label=f'R2 = {R2}')
            plt.legend()


def Colors(clusters):
    assign = []
    cols = ['r', 'b', 'g', 'c', 'orange', 'y', 'magenta' 'k']
    for j in clusters:
        assign.append(cols[j])
    return assign
```

9ee13e8e929b54ff11c55d633328677482df68a2
ipynb
Jupyter Notebook
Mate 2 Funciones.ipynb
lucas-abbate/mate2-bio
a0c70b9d0ef9bb183318bada8c5880c6c88ff1b5
MIT
Mate 2 Funciones.ipynb
lucas-abbate/mate2-bio
a0c70b9d0ef9bb183318bada8c5880c6c88ff1b5
MIT
Mate 2 Funciones.ipynb
lucas-abbate/mate2-bio
a0c70b9d0ef9bb183318bada8c5880c6c88ff1b5
MIT
Qwen/Qwen-72B
1. YES
2. YES