Pular para conteúdo

Lista 15 · Ajuste não linear

Acompanha o capítulo 15. Polinômios, linearização (exponencial e potência) e ajuste direto com curve_fit — e o cuidado com o sobreajuste.

Abrir no Colab

Baixar o caderno

Como esta lista funciona

Abra o caderno no Colab, escreva a sua solução no lugar do # TODO e rode a célula de teste logo abaixo. Ela mesma diz o que passou e o que falhou — você não precisa esperar correção de ninguém.

O código de cada exercício também está aqui embaixo, caso você prefira colar num editor local. Os testes comparam números com tolerância.

Exercícios marcados com ✏️ são à mão, no papel, com 4 casas decimais — como na parte em papel das avaliações.

Exercício 01

✏️ À mão. Uma colônia de bactérias foi contada de hora em hora: \((0, 5)\), \((1, 10)\), \((2, 20)\) (milhares).

a) Transforme os dados em \((x, \ln y)\).

b) Os três pontos transformados estão numa reta? Calcule a inclinação e o intercepto.

c) Escreva o modelo \(y = a\,e^{bx}\) e o tempo que a colônia leva para dobrar.

Exercício 02

Escreva ajuste_exponencial(x, y), que ajusta \(y = a\,e^{bx}\) pela linearização (\(\ln y\) contra \(x\), com np.polyfit) e devolve a tupla (a, b).

import numpy as np


def ajuste_exponencial(x, y):
    """(a, b) de y = a * e^(b x), pela reta de ln(y) contra x."""
    # TODO: np.polyfit no ln(y); lembre que o polyfit devolve [inclinação, intercepto]
    pass


# --- testes (não altere) ---
import numpy as np


r = ajuste_exponencial(np.array([0.0, 1.0, 2.0]), np.array([5.0, 10.0, 20.0]))
assert r is not None, "faltou o return?"
assert abs(r[0] - 5) < 1e-9 and abs(r[1] - np.log(2)) < 1e-9, f"veio {r}"
print("Exercício 02: todos os testes passaram!")

Exercício 03

Escreva ajuste_potencia(x, y), que ajusta \(y = a\,x^b\) pela linearização log-log e devolve (a, b).

import numpy as np


def ajuste_potencia(x, y):
    """(a, b) de y = a * x^b, pela reta de ln(y) contra ln(x)."""
    # TODO: np.polyfit com os dois logaritmos
    pass


# --- testes (não altere) ---
import numpy as np


r = ajuste_potencia(np.array([1.0, 2.0, 4.0]), np.array([3.0, 12.0, 48.0]))
assert r is not None, "faltou o return?"
assert abs(r[0] - 3) < 1e-9 and abs(r[1] - 2) < 1e-9, f"y = 3 x²: veio {r}"
print("Exercício 03: todos os testes passaram!")

Exercício 04

Medicina nuclear — a meia-vida. Escreva meia_vida(dias, contagens), que usa o ajuste_exponencial (já escrito) e devolve a meia-vida, \(\ln 2 / k\), com \(k = -b\).

meia_vida([0, 2, 4, 6, 8, 10, 12, 14, 16], [10020, 8440, 7010, 5950, 4970, 4230, 3500, 2980, 2510])    # -> ≈ 8.00

import numpy as np


def ajuste_exponencial(x, y):
    """JÁ ESCRITA."""
    coef = np.polyfit(x, np.log(y), 1)
    return (np.exp(coef[1]), coef[0])


def meia_vida(dias, contagens):
    """Meia-vida (dias) de uma amostra radioativa, pelo ajuste exponencial."""
    # TODO: o ajuste e a fórmula (cuidado com o sinal: b é negativo)
    pass


# --- testes (não altere) ---
import numpy as np


assert abs(meia_vida(np.array([0, 2, 4, 6, 8, 10, 12, 14, 16]), np.array([10020, 8440, 7010, 5950, 4970, 4230, 3500, 2980, 2510])) - 8.00499724700046) < 1e-9
assert abs(meia_vida(np.array([0.0, 5.0, 10.0]), np.array([80.0, 40.0, 20.0])) - 5) < 1e-9, "cai à metade a cada 5 dias"
print("Exercício 04: todos os testes passaram!")

Exercício 05

Escolhendo a complexidade. A função erro_rms(coef, x, y) (já escrita) mede o erro de um polinômio do polyfit. Escreva melhor_grau(x_treino, y_treino, x_teste, y_teste, graus), que ajusta um polinômio de cada grau da lista no treino e devolve o grau com o menor erro no teste.

import numpy as np


def erro_rms(coef, x, y):
    """Raiz do erro quadrático médio do polinômio coef (formato do polyfit). JÁ ESCRITA."""
    previsto = np.zeros(len(x))
    grau = len(coef) - 1
    for k in range(len(coef)):
        previsto = previsto + coef[k] * x**(grau - k)
    return np.sqrt(np.mean((previsto - y)**2))


def melhor_grau(x_treino, y_treino, x_teste, y_teste, graus):
    """O grau (da lista graus) com o MENOR erro no teste."""
    # TODO: padrão extremo sobre os graus; o erro é calculado nos dados de TESTE
    pass


# --- testes (não altere) ---
import numpy as np


x = np.linspace(0, 4, 21)
y = np.array([1.3, 1.18, 1.82, 2.27, 1.9800000000000002, 2.5, 2.8800000000000003, 2.72, 2.67, 3.1299999999999994, 3.05, 2.78, 3.2199999999999998, 2.77, 2.78, 2.2, 2.4799999999999995, 2.0199999999999996, 1.569999999999999, 1.6300000000000008, 0.9])
g = melhor_grau(x[0::2], y[0::2], x[1::2], y[1::2], [1, 2, 3, 5, 7, 9])
assert g == 2, f"veio {g}"
print("Exercício 05: todos os testes passaram!")

Exercício 06

Indústria — a curva de aprendizagem. Numa linha de montagem (a Boeing mediu isso nos anos 1930), cada vez que a produção dobra, o tempo para montar uma unidade cai para uma fração fixa: \(T(n) = T_1\,n^b\), com \(b < 0\). Escreva taxa_de_aprendizagem(unidades, tempos), que ajusta a potência (o ajuste_potencia já está escrito) e devolve \(2^b\).

taxa_de_aprendizagem([1, 2, 4, 8, 16, 32], [120.0, 98.4, 80.7, 66.2, 54.3, 44.5])    # -> ≈ 0.820

Isso quer dizer uma curva de aprendizagem de quantos por cento?

import numpy as np


def ajuste_potencia(x, y):
    """JÁ ESCRITA."""
    coef = np.polyfit(np.log(x), np.log(y), 1)
    return (np.exp(coef[1]), coef[0])


def taxa_de_aprendizagem(unidades, tempos):
    """Fração do tempo que sobra quando a produção dobra: 2 elevado a b."""
    # TODO: o ajuste e 2**b
    pass


# --- testes (não altere) ---
import numpy as np


assert abs(taxa_de_aprendizagem(np.array([1, 2, 4, 8, 16, 32]), np.array([120.0, 98.4, 80.7, 66.2, 54.3, 44.5])) - 0.8200975212199776) < 1e-9
assert abs(taxa_de_aprendizagem(np.array([1.0, 2.0, 4.0]), np.array([100.0, 90.0, 81.0])) - 0.9) < 1e-9
print("Exercício 06: todos os testes passaram!")

Exercício 07

Saúde — o tamanho de uma epidemia. Escreva total_da_epidemia(t, casos), que ajusta a logistica (já escrita) com curve_fit e devolve o \(K\), o total de casos quando a epidemia acabar. Use como chutes iniciais p0=[np.max(casos), 50, 0.1].

Com dados só até o pico da epidemia, o \(K\) ajustado é confiável? Teste cortando os dados pela metade.

import numpy as np
from scipy.optimize import curve_fit


def logistica(t, K, A, r):
    """JÁ ESCRITA."""
    return K / (1 + A * np.exp(-r * t))


def total_da_epidemia(t, casos):
    """O K (total de casos) da logística ajustada com curve_fit."""
    # TODO: curve_fit(logistica, t, casos, p0=[...]) e o primeiro parâmetro
    pass


# --- testes (não altere) ---
import numpy as np


t = np.linspace(0, 60, 31)
casos = logistica(t, 5000.0, 80.0, 0.15)
assert abs(total_da_epidemia(t, casos) - 5000) < 1e-3, "dados exatos: recupera o K"
print("Exercício 07: todos os testes passaram!")