How to modify the value after a judge

Questions
1

Now I have a code as follows


import matplotlib.pyplot as plt
import numpy as np
from scipy import stats
import seaborn as sns
import arviz as az
import numpy as np
import pymc3 as pm
import os
import theano.tensor as tt
from theano.compile.ops import as_op
import theano

no_data=100

x=np.random.rand(no_data)
y=np.random.rand(no_data)


a=5
b=10
z=a*x+b*y

for i in range(no_data):
    if z[i]<4:
        z[i]=0



test_model = pm.Model()


with test_model:
    a_fit = pm.Categorical('a_fit', p=np.ones(10)/10)
    b_fit = pm.Categorical('b_fit', p=np.ones(20)/20)  
    
    z_fit=a_fit*x+b_fit*y
    for j in range(no_data):
        if tt.lt(z_fit[j],4):  #z_fit[j]<4:
            z_fit=tt.set_subtensor(z_fit[j],0)    #z_fit[j]=0   
    
    z_obs_nor=pm.Normal("z_fit", mu=z_fit, sigma=1, observed=z)
    trace_nor=pm.sample(draws=10000,tune=100,chains=4,cores=4)
    
pm.traceplot(trace_nor)

This code give a wrong estimation of a and b.
Can anyone give me some suggestion?
thanks a lot.

2

Hi, your issue is that you are using Python’s standard if statement rather than using one of Theano’s logical operators. See this collab notebook for a detailed explanation, but here’s my suggested fix (avoiding both an if and your loop altogether):

test_model = pm.Model()

with test_model:
    a_fit = pm.Categorical('a_fit', p=np.ones(10)/10)
    b_fit = pm.Categorical('b_fit', p=np.ones(20)/20)  
    
    z_fit = a_fit*x+b_fit*y
    
    # Get index values where entry is < 4
    idx_pts = (z_fit < 4).nonzero()

    # Set these entries from above to 0
    z_fit = tt.set_subtensor(z_fit[idx_pts], 0)
    
    z_obs_nor = pm.Normal("z_fit", mu=z_fit, sigma=1, observed=z)
    trace_nor=pm.sample(draws=4000,tune=1000,chains=4,cores=4)

This correctly estimates a and b. Hope this helps!

1 Like
3

It worked.
Thanks for your reply.

1 Like