Priesemann-Group / covid19_inference

Bayesian python toolbox for inference and forecast of the spread of the Coronavirus
GNU General Public License v3.0
73 stars 70 forks source link

Example SIR model with one german state has negative lambda values? #32

Closed anshu0612 closed 4 years ago

anshu0612 commented 4 years ago

Why the values for the lambda_0, lambda_1, lambda_2, and lambda_3 negative. Aren't they supposed to be positive, at least in the initial stages when the number of cases is indeed increasing? I tried for Singapore and Itay and got similar results.

What am I missing here?

semohr commented 4 years ago

I just tried to reproduce this with the jhu dataset for Italy over the period from 2020/3/14 to today. A quick run with 500 draws seem to work fine for me. Would you maybe give us some code snippets of the main things you changed in the example?

Additional what are your changepoint priors?

I get the following graph for lambda_t without changing anything except the data retrieval in the example.

lambda_t_italy

y,x=cov19.plot._get_array_from_trace_via_date(model,trace,"lambda_t",bd,model.sim_end) ax=cov19.plot._timeseries(x,y,what="fcast",draw_ci_95=True)
cov19.plot._format_date_xticks(ax) ax.set_ylabel(r"$\lambda$")

semohr commented 4 years ago

If you are looking at the lambda_0_log, lambda_1_log, lambda_2_log distributions it is very much possible that they are negative. For example log(0.1)=-1.

anshu0612 commented 4 years ago

hi @semohr

I tried running the model from 4th March to 6th May on JHU dataset. I took three policy priors:

# All sporting event canceled in Italy
prior_date_mild_dist_begin =  datetime.datetime(2020,3,9)
# Nation wide lockdown
prior_date_strong_dist_begin =  datetime.datetime(2020,3,21)
#Lockdown extended tiil  3 May
prior_date_contact_ban_begin =  datetime.datetime(2020,4,10)

Here is the link to the notebook.
https://github.com/anshu0612/covid19_inference/blob/master/scripts/Italy.ipynb

So I need to take the log of the x-axis values in order to get the Lambda?

The documentation says that the notebook example provided on the single state is similar to the paper, but I can't see any similarity

semohr commented 4 years ago

If my understanding is correct, taking the log on the x-axis should give you your expected value. Maybe you want to look at the "lambda_t" trace as well.

In your code ax.set_xlabel(varnames[i_ax].replace("_log", "")) replaces the log string. Maybe that was reason for the confusion.

I'm not too sure but the model should be the same as the one used in the paper. Maybe @jdehning can give more insight on the changes in the model (if there are any).

anshu0612 commented 4 years ago

@semohr Yeah I added ax.set_xlabel(varnames[i_ax].replace("_log", "")). That's not the problem I can't take log along the x-axis, taking log of the negative values will result in complex number. And yes I checked the traces, they all are negative values

Konrad982 commented 4 years ago

So if I understand correctly, you plot the histograms of the logarithms of the lambdas and then you exchange the label to show only "lambda_i". I think you need to take the exponential of the ticks of the x-axis.

You may also plot directly a histogram of the lambdas modifying plotting.py, as I did here. I did however cut off the data above and below some quantiles in order maintain readability.

anshu0612 commented 4 years ago

Hi @Konrad982 My concern is to get rid of the negative values, or get an understanding of the negative values, why all lambdas are negative and how to replace with the actual value which ranges from 0 to 1.

Hope so yo get my concern.

Konrad982 commented 4 years ago

What you see there should be the histogram of the logarithm of lambda (which indeed may be negative), even if the label tells you that it is lambda. That's why you need to take the exponential of the values on the x-axis in order to get rid of the negative values.

anshu0612 commented 4 years ago

@Konrad982 Thanks :)