Hello, Thanks for your effort. I'm verifying your code.

In model.test() in net.py, decoding part seems unnatural.

The task is to predict indices from 168 ~ 191 (24 time samples) using 0 ~ 167 samples (168 time samples). It is assumed that we do not know predicting part.

In your source, however, at iter 0, we already know 168 index from original data, and it occurred last time step to be discarded due to out of index.

I thinks it can be fixed like this:

for t in range(self.params.predict_steps):
   mu_de, sigma_de, decoder_hidden, decoder_cell = self(x[self.params.predict_start+(t-1)].unsqueeze(0), id_batch, decoder_hidden, decoder_cell)
   sample_mu[:, t] = mu_de * v_batch[:, 0] + v_batch[:, 1]
   sample_sigma[:, t] = sigma_de * v_batch[:, 0]
   x[self.params.predict_start+t, :, 0] = mu_de

Hello, thank you that i have implemented the DeepAR model for my data set. but i want to implement negative Binomial distributions likelihood loss function since it is in the paper. I am very beginner in PyTorch but somehow I have implemented your Gaussian implementation in for my dataset. Pytorch itself has a Negative Binomial Class which contains total count, logits, probs as a parameter. actually I cant relate which of these parameters indicates the zero index, mu, sigma I have tried with my intuitions. but my y_prediction output is showing inappropriate value. I am looking for your suggestion which might be helpful for me to solve this

In the /model/net.py file, line 87, when the model do inference from [predict_start ,predict_start+ predict_step] , it can use the true information of them, which is not properly for model test. In the test stage, the model only input the true value of x[predict_start], and the following step should input the last predict mu instead.

Hello. Thanks so much for creating this repository. I was going through the preprocessing code, and it seems like there is an error in the code. But I wanted to check with you, in case I was just misinterpreting it.

So in the preprocess_elect.py module, in the prep_data() function there is a set of lines like:

            x_input[count, 1:, 0] = data[window_start:window_end-1, series]  # <--- QUESTION
            x_input[count, :, 1:1+num_covariates] = covariates[window_start:window_end, :]
            x_input[count, :, -1] = series

The problem is that the x_input[count, 1:, 0] begins at index 1 instead of 0. So that will actually make the window size 191 instead of 192. And it is not clear what the benefit of starting every window with 0 will be. But further, the other arrays, such as labels goes from window_start to window_end. If that is the case, then the x_input array and the labels array is off by an increment.

So I was not sure which is the right way to do it, or whether this is written accurately. Would you mind clarifying if you have a chance. Thanks.

Hi there, thanks for sharing the code. Can you also elaborate on the results/graphs that you have plotted. I believe the blue line is the actual time series data(ground truth) and the red line is our Gaussian likelihood function for which we are also measuring those metrics like ND and RMSE. This is done on 200 samples/days from the time series at a time, is my understanding correct? In my case, I have taken 150 days of data(130 days train and 20 days test) , so I don't understand how the graph I got after training has portions after 150, will you please expain. Thanks!

first，thanks for ur code！it‘s very helpful! in your code， when your calculate the mu&sigma use the h from all three layers。

but in the deepar paper calculate the mu&sigma use the h from the last layers 。

which one is better？

the signature of gen_covariates in preprocess_elect.py takes a parameter num_covariates, and uses it to reshape the array it returns. However, this will crash for num_covariates < 4, since the shape of the containing array is set by num_covariates and column indices 1, 2, 3 are accessed while it is being populated with data. If the containing array is (n, 3), the function will crash when it tries to put the month in covariates[i, 3].

def gen_covariates(times, num_covariates):
    covariates = np.zeros((times.shape[0], num_covariates))
    for i, input_time in enumerate(times):
        covariates[i, 1] = input_time.weekday()
        covariates[i, 2] = input_time.hour
        covariates[i, 3] = input_time.month
    for i in range(1,num_covariates):
        covariates[:,i] = stats.zscore(covariates[:,i])
    return covariates[:, :num_covariates]

Hello, I am new in PyTorch but i implemented DeepAR model by your implementation. I need a help to make time window inference. i am saving model as .pkl in the base model and I have modified the evaluate.py code a bit. what i need is predict function that will take a time window and give the prediction result. so far I have done this

# -*- coding: utf-8 -*-
"""
Created on Sat Jun 13 13:42:11 2020

"""

import argparse
import logging
import os

import numpy as np
import torch
from torch.utils.data.sampler import RandomSampler
from tqdm import tqdm

import utils
import model.net as net
from dataloader import *

import matplotlib
matplotlib.use('Agg')
import matplotlib.pyplot as plt
import pandas as pd

logger = logging.getLogger('DeepAR.Eval')

parser = argparse.ArgumentParser()
parser.add_argument('--dataset', default='test_dataset', help='Name of the dataset')
parser.add_argument('--data-folder', default='data', help='Parent dir of the dataset')
parser.add_argument('--model-name', default='base_model', help='Directory containing params.json')
parser.add_argument('--relative-metrics', action='store_true', help='Whether to normalize the metrics by label scales')
parser.add_argument('--sampling', action='store_true', help='Whether to sample during evaluation')
parser.add_argument('--restore-file', default='best',
                    help='Optional, name of the file in --model_dir containing weights to reload before \
                    training')



def evaluate(model, loss_fn, test_loader, params, plot_num, sample=True):
    '''Evaluate the model on the test set.
    Args:
        model: (torch.nn.Module) the Deep AR model
        loss_fn: a function that takes outputs and labels per timestep, and then computes the loss for the batch
        test_loader: load test data and labels
        params: (Params) hyperparameters
        plot_num: (-1): evaluation from evaluate.py; else (epoch): evaluation on epoch
        sample: (boolean) do ancestral sampling or directly use output mu from last time step
    '''
    model.eval()
    with torch.no_grad():
      print("In model .eval")  
      plot_batch = np.random.randint(len(test_loader))
      print("plot_batch value in evaluation", plot_batch)

      summary_metric = {}
      raw_metrics = utils.init_metrics(sample=sample)

      # Test_loader: 
      # test_batch ([batch_size, train_window, 1+cov_dim]): z_{0:T-1} + x_{1:T}, note that z_0 = 0;
      # id_batch ([batch_size]): one integer denoting the time series id;
      # v ([batch_size, 2]): scaling factor for each window;
      # labels ([batch_size, train_window]): z_{1:T}.
      for i, (test_batch, id_batch, v, labels) in enumerate(tqdm(test_loader)):
          test_batch = test_batch.permute(1, 0, 2).to(torch.float32).to(params.device)
          id_batch = id_batch.unsqueeze(0).to(params.device)
          v_batch = v.to(torch.float32).to(params.device)
          labels = labels.to(torch.float32).to(params.device)
          batch_size = test_batch.shape[1]
          input_mu = torch.zeros(batch_size, params.test_predict_start, device=params.device) # scaled
          input_sigma = torch.zeros(batch_size, params.test_predict_start, device=params.device) # scaled
          hidden = model.init_hidden(batch_size)
          cell = model.init_cell(batch_size)

          for t in range(params.test_predict_start):
              # if z_t is missing, replace it by output mu from the last time step
              zero_index = (test_batch[t,:,0] == 0)
              if t > 0 and torch.sum(zero_index) > 0:
                  test_batch[t,zero_index,0] = mu[zero_index]

              mu, sigma, hidden, cell = model(test_batch[t].unsqueeze(0), id_batch, hidden, cell)
              input_mu[:,t] = v_batch[:, 0] * mu + v_batch[:, 1]
              input_sigma[:,t] = v_batch[:, 0] * sigma

          if sample:
              samples, sample_mu, sample_sigma = model.test(test_batch, v_batch, id_batch, hidden, cell, sampling=True)
              raw_metrics = utils.update_metrics(raw_metrics, input_mu, input_sigma, sample_mu, labels, params.test_predict_start, samples, relative = params.relative_metrics)
          else:
              sample_mu, sample_sigma = model.test(test_batch, v_batch, id_batch, hidden, cell)
              raw_metrics = utils.update_metrics(raw_metrics, input_mu, input_sigma, sample_mu, labels, params.test_predict_start, relative = params.relative_metrics)

          if i == plot_batch:
              if sample:
                  sample_metrics = utils.get_metrics(sample_mu, labels, params.test_predict_start, samples, relative = params.relative_metrics)
              else:
                  sample_metrics = utils.get_metrics(sample_mu, labels, params.test_predict_start, relative = params.relative_metrics)                
              # select 10 from samples with highest error and 10 from the rest
              top_10_nd_sample = (-sample_metrics['ND']).argsort()[:batch_size // 10]  # hard coded to be 10
              chosen = set(top_10_nd_sample.tolist())
              all_samples = set(range(batch_size))
              not_chosen = np.asarray(list(all_samples - chosen))
              if batch_size < 100: # make sure there are enough unique samples to choose top 10 from
                  random_sample_10 = np.random.choice(top_10_nd_sample, size=10, replace=True)
              else:
                  random_sample_10 = np.random.choice(top_10_nd_sample, size=10, replace=False)
              if batch_size < 12: # make sure there are enough unique samples to choose bottom 90 from
                  random_sample_90 = np.random.choice(not_chosen, size=10, replace=True)
              else:
                  random_sample_90 = np.random.choice(not_chosen, size=10, replace=False)
              combined_sample = np.concatenate((random_sample_10, random_sample_90))

              label_plot = labels[combined_sample].data.cpu().numpy()
              predict_mu = sample_mu[combined_sample].data.cpu().numpy()
              predict_sigma = sample_sigma[combined_sample].data.cpu().numpy()
              plot_mu = np.concatenate((input_mu[combined_sample].data.cpu().numpy(), predict_mu), axis=1)
              plot_sigma = np.concatenate((input_sigma[combined_sample].data.cpu().numpy(), predict_sigma), axis=1)
              plot_metrics = {_k: _v[combined_sample] for _k, _v in sample_metrics.items()}
              plot_eight_windows(params.plot_dir, plot_mu, plot_sigma, label_plot, params.test_window, params.test_predict_start, plot_num, plot_metrics, sample)

      summary_metric = utils.final_metrics(raw_metrics, sampling=sample)
      metrics_string = '; '.join('{}: {:05.3f}'.format(k, v) for k, v in summary_metric.items())
      logger.info('- Full test metrics: ' + metrics_string)
    return summary_metric




def plot_eight_windows(plot_dir,
                       predict_values,
                       predict_sigma,
                       labels,
                       window_size,
                       predict_start,
                       plot_num,
                       plot_metrics,
                       sampling=False):

    x = np.arange(window_size)
    f = plt.figure(figsize=(8, 42), constrained_layout=True)
    nrows = 21
    ncols = 1
    ax = f.subplots(nrows, ncols)

    for k in range(nrows):
        if k == 10:
            ax[k].plot(x, x, color='g')
            ax[k].plot(x, x[::-1], color='g')
            ax[k].set_title('This separates top 10 and bottom 90', fontsize=10)
            continue
        m = k if k < 10 else k - 1
        # for 24 hour predictioon
        
        print('iteration',k)
        predict_or_df=pd.DataFrame(predict_values[m,predict_start:],columns=['y_pred'])
        predict_next_df=pd.DataFrame(labels[m,predict_start:],columns=['y_true'])
        y_lower= predict_values[m,predict_start:]-predict_sigma[m,predict_start:]
        y_lower_df=pd.DataFrame(y_lower,columns=['y_lower'])
        y_upper=predict_values[m,predict_start:]+ predict_sigma[m,predict_start:]
        y_upper_df=pd.DataFrame(y_upper,columns=['y_upper'])

        y2_lower=predict_values[m,predict_start:] - 2 * predict_sigma[m,predict_start:]
        y2_lower_df=pd.DataFrame(y2_lower,columns=['y2_lower'])
        y2_upper=predict_values[m,predict_start:] + 2 * predict_sigma[m,predict_start:]
        y2_upper_df=pd.DataFrame(y2_upper,columns=['y2_upper'])

        frames=[predict_next_df,predict_or_df,y_lower_df,y_upper_df,y2_lower_df,y2_upper_df]
        summarydf=pd.concat(frames,axis=1)
        print(summarydf)
        summarydf.to_csv('/Users/kalyan.admin/KalyanWork/Timeseries_covid/experiments/csvfiles_inferenceset/pytorch_report_"%s".csv'%(k))















       #Visualgraph for 192 hour windows that is  week based hourly prediction
        ax[k].plot(x, predict_values[m], color='b')
        ax[k].fill_between(x[predict_start:], predict_values[m, predict_start:] - 2 * predict_sigma[m, predict_start:],
                         predict_values[m, predict_start:] + 2 * predict_sigma[m, predict_start:], color='blue',
                         alpha=0.2)
        ax[k].plot(x, labels[m, :], color='r')
        ax[k].axvline(predict_start, color='g', linestyle='dashed')

        #metrics = utils.final_metrics_({_k: [_i[k] for _i in _v] for _k, _v in plot_metrics.items()})


        plot_metrics_str = f'ND: {plot_metrics["ND"][m]: .3f} ' \
            f'RMSE: {plot_metrics["RMSE"][m]: .3f}'
        if sampling:
            plot_metrics_str += f' rou90: {plot_metrics["rou90"][m]: .3f} ' \
                                f'rou50: {plot_metrics["rou50"][m]: .3f}'

        ax[k].set_title(plot_metrics_str, fontsize=10)

    f.savefig(os.path.join(plot_dir, str(plot_num) + '.png'))
    plt.close()



if __name__ == '__main__':
    # Load the parameters
    args = parser.parse_args()
    model_dir = os.path.join('experiments', args.model_name) 
    json_path = os.path.join(model_dir, 'params.json')
    data_dir = os.path.join(args.data_folder, args.dataset)
    assert os.path.isfile(json_path), 'No json configuration file found at {}'.format(json_path)
    params = utils.Params(json_path)

    utils.set_logger(os.path.join(model_dir, 'eval.log'))

    params.relative_metrics = args.relative_metrics
    params.sampling = args.sampling
    params.model_dir = model_dir
    params.plot_dir = os.path.join(model_dir, 'figures')
    
    cuda_exist = torch.cuda.is_available()  # use GPU is available

    # Set random seeds for reproducible experiments if necessary
    if cuda_exist:
        params.device = torch.device('cuda')
        # torch.cuda.manual_seed(240)
        logger.info('Using Cuda...')
        model = net.Net(params).cuda()
    else:
        params.device = torch.device('cpu')
        # torch.manual_seed(230)
        logger.info('Not using cuda...')
        model = net.Net(params)

    # Create the input data pipeline
    logger.info('Loading the datasets...')

    test_set = TestDataset(data_dir, args.dataset, params.num_class)
    test_loader = DataLoader(test_set, batch_size=params.predict_batch, sampler=RandomSampler(test_set), num_workers=4)
    print("testloader",type(test_loader))
    #print(test_loader)
    logger.info('- done.')

    print('model: ', model)
    loss_fn = net.loss_fn

    logger.info('Starting evaluation')

    # Reload weights from the saved file
    utils.load_checkpoint(os.path.join(model_dir, args.restore_file + '.pkl'), model)

    test_metrics = evaluate(model, loss_fn, test_loader, params, -1, params.sampling)
    save_path = os.path.join(model_dir, 'metrics_test_{}.json'.format(args.restore_file))
    utils.save_dict_to_json(test_metrics, save_path)

so up to this in plot window function plot_eight_windows i am calculating my y_pred, y_true, y_upper,y_lower, and save those files into a folder. what i have done is i have processed a test dataset and saved it in differnt folder for only inference purposes. and i use this it is in pretrained model in best_pkl. I am getting result but i am in confusion like , am i evaluating here or doing actual inference. what if I call a time range

2020-04-25 00:00:00 and want to predict the related values o. using my pretrained model. Kindly help me to sort it out

@zhykoties Hi there, I just wanted to understand the uncertainty estimation graph. Here, red line denotes the likelihood function and the blue denotes the actual ground truth data. In that case, the X-axis is the sample size of the data which is 200 here but I did not understand the Y-axis, what does the Y-axis represent exactly, will you please explain. Also, the values change in the Y-axis. Thanks!

hello, why is it that v_input is initialized with 2 columns (mean and std i presume) but only one of them is changed. It also looks strange that instead of subtracting the mean, the series is divided by the mean as in the following lines:

v_input[count, 0] = np.true_divide(x_input[count, 1:input_size, 0].sum(),nonzero_sum)+1
x_input[count, :, 0] = x_input[count, :, 0]/v_input[count, 0]

I would like to access to the preprocessing realised on these others dataset that appears on the original papers to validate the results, is it possible? Regards.

for different train test split dataset I am getting this error

Traceback (most recent call last):

  File "<ipython-input-4-6b0715d9506b>", line 6, in <module>
    args.restore_file)

  File "<ipython-input-1-a393da520de1>", line 121, in train_and_evaluate
    test_loader, params, epoch)

  File "<ipython-input-1-a393da520de1>", line 83, in train
    test_metrics = evaluate(model, loss_fn, test_loader, params, epoch, sample=args.sampling)

  File "/Users/kalyan.admin/KalyanWork/Timeseries_covid/evaluate.py", line 95, in evaluate
    random_sample_10 = np.random.choice(top_10_nd_sample, size=10, replace=True)

  File "mtrand.pyx", line 773, in numpy.random.mtrand.RandomState.choice

the implementation worked well. But right now I am getting it for several train test splitting. what will be the solution