Loss constant

Hi !

I’m trying to work on MLE and MAP with pyro for my internship.
Here is my issue :
I’ve got 3 parameter (I0, beta, T) that I would like to maximize.
However my program doesn’t run as it should the optimizer : loss doesn’t reduce over iterations.

import math
import numpy as np
import matplotlib.pyplot as plt
import os
import torch
import torch.distributions.constraints as constraints
import pyro
from pyro.optim import Adam
from pyro.infer import SVI, Trace_ELBO
import pyro.distributions as dist

this is for running the notebook in our testing framework

smoke_test = “CI” in os.environ
n_steps = 2 if smoke_test else 2000

assert pyro.version.startswith(“1.8.1”)

clear the param store in case we’re in a REPL

pyro.clear_param_store()

Lambda_250 = torch.tensor(250e-6)
Lambda_360 = torch.tensor(360e-6)
Lambda_520 = torch.tensor(520e-6)

I0_vrai = torch.tensor(1.0)
beta_vrai = torch.tensor(1.0)
T_vrai = torch.tensor(18.0)

def f_beta(Lambda, T):

k = np.float64(1.38e-23)  # constante de Boltztmann
h = np.float64(6.626e-34)  # constante de Planck
c = np.float64(3e8)  # célérité

return (
    2
    * h
    * pow(c, 2)
    / pow(Lambda, 5)
    * 1
    / (torch.exp(h * c / (Lambda * k * T)) - 1)
)

def I(I0, beta, T, Lambda):
return I0 * torch.pow(Lambda, -beta) * f_beta(Lambda, T)

def H(I0, beta, T):
return torch.tensor(
[
I(I0, beta, T, Lambda_250),
I(I0, beta, T, Lambda_360),
I(I0, beta, T, Lambda_520),
]
)

data = torch.tensor(
[
I(I0_vrai, beta_vrai, T_vrai, Lambda_250),
I(I0_vrai, beta_vrai, T_vrai, Lambda_360),
I(I0_vrai, beta_vrai, T_vrai, Lambda_520),
]
)

noise = torch.div(data, 100, rounding_mode=“trunc”)

data bruitées

data = torch.tensor(
[
dist.Normal(I(I0_vrai, beta_vrai, T_vrai, Lambda_250), noise[0]).sample(),
dist.Normal(I(I0_vrai, beta_vrai, T_vrai, Lambda_360), noise[1]).sample(),
dist.Normal(I(I0_vrai, beta_vrai, T_vrai, Lambda_520), noise[2]).sample(),
]
)

écart type des lois normales

std = torch.tensor(100.0)

Déinition du modèle

def model(data):

# sample f from the Uniform prior
I0 = pyro.sample("I0", dist.Uniform(0, 5))
beta = pyro.sample("beta", dist.Uniform(0, 2))
T = pyro.sample("T", dist.Uniform(5, 100))

# loop over the observed data
with pyro.plate("data", data.size(0)):
    pyro.sample("obs", dist.Normal(H(I0, beta, T), std), obs=data)

def guide(data):

# register the three variational parameters with Pyro.

I0_q = pyro.param(“I0_q”, torch.tensor(1), constraint=constraints.positive)

beta_q = pyro.param(

“beta_q”, torch.tensor(1), constraint=constraints.interval(0.0, 2.0)

)

T_q = pyro.param(

“T_q”, torch.tensor(18), constraint=constraints.interval(5, 100)

)

# sample latent_intesity from the distribution Normale(I(I0_q, beta_q, T_q)

pyro.sample(“latent_intensity”, dist.Normal(H(I0_q, beta_q, T_q), std))

autoguide = pyro.infer.autoguide.AutoDelta(model)

Set up de la stochastic variational inference (svi)

Set up de la stochastic variational inference (svi)

def train(model, guide, lr=0.005, n_steps=201):
pyro.clear_param_store()
adam_params = {“lr”: lr}
adam = pyro.optim.Adam(adam_params)
svi = SVI(model, guide, adam, loss=Trace_ELBO())

for step in range(n_steps):
    loss = svi.step(data)
    if step % 50 == 0:
        print("[iter {}]  loss: {:.4f}".format(step, loss))

train(model, autoguide)

grab the learned variational parameters

I0_q = pyro.param(“I0_q”).item()

beta_q = pyro.param(“beta_q”).item()

T_q = pyro.param(“T_q”).item()

print(I0_q, beta_q, T_q)

print(
“Our estimated parameter I0 is {:.3f}”.format(autoguide.median(data)[“I0”].item())
)
print(
“Our estimated parameter beta is {:.3f}”.format(
autoguide.median(data)[“beta”].item()
)
)
print(“Our estimated parameter T is {:.3f}”.format(autoguide.median(data)[“T”].item()))

I0_est = torch.tensor(autoguide.median(data)[“I0”].item())
beta_est = torch.tensor(autoguide.median(data)[“beta”].item())
T_est = torch.tensor(autoguide.median(data)[“T”].item())

def model_mle(data):

# sample f from the Uniform prior
# I0_mle = pyro.param("I0_mle", I0_est, constraint=constraints.positive,)
# beta_mle = pyro.param("beta_mle", beta_est, constraint=constraints.positive,)
# T_mle = pyro.param("T_mle", T_est, constraint=constraints.positive,)

I0_mle = pyro.param("I0_mle", torch.tensor(2.0), constraint=constraints.positive,)
beta_mle = pyro.param(
    "beta_mle", torch.tensor(0.8), constraint=constraints.positive,
)
T_mle = pyro.param("T_mle", torch.tensor(35.0), constraint=constraints.positive,)

# loop over the observed data
with pyro.plate("data", data.size(0)):
    pyro.sample("obs", dist.Normal(H(I0_mle, beta_mle, T_mle), std), obs=data)

def guide_mle(data):
pass

train(model_mle, guide_mle)

mle_estimate_I0 = pyro.param(“I0_mle”).item()
mle_estimate_beta = pyro.param(“beta_mle”).item()
mle_estimate_T = pyro.param(“T_mle”).item()
print(“Our MLE estimate of the I0 parameter is {:.3f}”.format(mle_estimate_I0))
print(“Our MLE estimate of the beta parameter is {:.3f}”.format(mle_estimate_beta))
print(“Our MLE estimate of the T parameter is {:.3f}”.format(mle_estimate_T))

Have you any idea ?

hello your code is hard to read because of lack of formatting (backticks ``` etc) so i haven’t read your code in detail but i suggest you read this tutorial, paying particular attention to sections 7 and 11

I’m trying to insert it more clearly :

My function that i use just below :

pyro.clear_param_store()

Lambda_250 = torch.tensor(250e-6)
Lambda_360 = torch.tensor(360e-6)
Lambda_520 = torch.tensor(520e-6)

I0_vrai = torch.tensor(1.0)
beta_vrai = torch.tensor(1.0)
T_vrai = torch.tensor(18.0)

def f_beta(Lambda, T):

k = np.float64(1.38e-23)  # constante de Boltztmann
h = np.float64(6.626e-34)  # constante de Planck
c = np.float64(3e8)  # célérité

return (
    2
    * h
    * pow(c, 2)
    / pow(Lambda, 5)
    * 1
    / (torch.exp(h * c / (Lambda * k * T)) - 1)
)

def I(I0, beta, T, Lambda):
return I0 * torch.pow(Lambda, -beta) * f_beta(Lambda, T)

def H(I0, beta, T):
return torch.tensor(
[
I(I0, beta, T, Lambda_250),
I(I0, beta, T, Lambda_360),
I(I0, beta, T, Lambda_520),
]
)

#Data noised used
data = torch.tensor(
[
I(I0_vrai, beta_vrai, T_vrai, Lambda_250),
I(I0_vrai, beta_vrai, T_vrai, Lambda_360),
I(I0_vrai, beta_vrai, T_vrai, Lambda_520),
]
)

noise = torch.div(data, 1000, rounding_mode=“trunc”)

data = torch.tensor(
[
dist.Normal(I(I0_vrai, beta_vrai, T_vrai, Lambda_250), noise[0]).sample(),
dist.Normal(I(I0_vrai, beta_vrai, T_vrai, Lambda_360), noise[1]).sample(),
dist.Normal(I(I0_vrai, beta_vrai, T_vrai, Lambda_520), noise[2]).sample(),
]
)

std = torch.tensor(0.05)

def model(data):

I0 = pyro.sample("I0", dist.Uniform(0, 5))
beta = pyro.sample("beta", dist.Uniform(0, 2))
T = pyro.sample("T", dist.Uniform(5, 100))

with pyro.plate("data", data.size(0)):
    pyro.sample("obs", dist.Normal(H(I0, beta, T), std), obs=data)

autoguide = pyro.infer.autoguide.AutoDelta(model)

def train(model, guide, lr=0.005, n_steps=201):
pyro.clear_param_store()
adam_params = {“lr”: lr}
adam = pyro.optim.Adam(adam_params)
svi = SVI(model, guide, adam, loss=Trace_ELBO())

for step in range(n_steps):
    loss = svi.step(data)
    if step % 50 == 0:
        print("[iter {}]  loss: {:.4f}".format(step, loss))

train(model, autoguide)

print(
“Our estimated parameter I0 is {:.3f}”.format(autoguide.median(data)[“I0”].item())
)
print(
“Our estimated parameter beta is {:.3f}”.format(
autoguide.median(data)[“beta”].item()
)
)
print(“Our estimated parameter T is {:.3f}”.format(autoguide.median(data)[“T”].item()))

Is it easier to read ?