Prev Next user_trace_init.py

@(@\newcommand{\B}[1]{ {\bf #1} } \newcommand{\R}[1]{ {\rm #1} } \newcommand{\W}[1]{ \; #1 \; }@)@This is dismod_at-20221105 documentation: Here is a link to its current documentation .
Using Initialization Trace Option

Purpose
Integrands
Nodes
True Iota
Model
     Parent Node
     Child Nodes
Data
trace_init_fit_model
Source Code

Purpose
This example shows how to use the trace_init_fit_model option.

Integrands
For this example there are one integrand, Sincidence.

Nodes
There three nodes. The first is called world and is the parent node for this example. The second (third) is called child_1 (child_2) and is a child of the parent node.

True Iota
For this example, the true model incidence rate iota is def iota_true(age, node) : import math iota_parent = 0.01 * (1 + age / 100.0) child_effect = 0.2 if node == 'world' : return iota_parent if node == 'child_1' : return math.exp(+ child_effect) * iota_parent if node == 'child_2' : return math.exp(- child_effect) * iota_parent assert False

Model

Parent Node
There is only one rate iota and it linear in age with knots at the age 0 and 100.

Child Nodes
There is only one rate iota and it constant.

Data
There are six data points measuring Sincidence with the true value of iota . These correspond to ages 0 and 100 at each of the three nodes.

trace_init_fit_model
The option trace_init_fit_model is set to true and the corresponding output is checked.

Source Code

# ------------------------------------------------------------------------
import sys
import os
import csv
import copy
import math
test_program = 'example/user/trace_init.py'
if sys.argv[0] != test_program  or len(sys.argv) != 1 :
   usage  = 'python3 ' + test_program + '\n'
   usage += 'where python3 is the python 3 program on your system\n'
   usage += 'and working directory is the dismod_at distribution directory\n'
   sys.exit(usage)
print(test_program)
#
# import dismod_at
local_dir = os.getcwd() + '/python'
if( os.path.isdir( local_dir + '/dismod_at' ) ) :
   sys.path.insert(0, local_dir)
import dismod_at
#
# change into the build/example/user directory
if not os.path.exists('build/example/user') :
   os.makedirs('build/example/user')
os.chdir('build/example/user')
# ------------------------------------------------------------------------
# Note that the a, t values are not used for this example
def example_db (file_name) :
   def fun_iota_parent(a, t) :
      return ('prior_iota_parent', None, None)
   def fun_iota_child(a, t) :
      return ('prior_iota_child', None, None)
   # ----------------------------------------------------------------------
   # age table:
   age_list    = [ 0.0, 100.0 ]
   #
   # time table:
   time_list   = [ 1990.0, 2020.0 ]
   #
   # integrand table:
   integrand_table = [
       { 'name':'Sincidence' }
   ]
   #
   # node table:
   node_table = [
      { 'name':'world',   'parent':''      } ,
      { 'name':'child_1', 'parent':'world' } ,
      { 'name':'child_2', 'parent':'world' } ,
   ]
   #
   # weight table:
   weight_table = list()
   #
   # covariate table:
   covariate_table = list()
   #
   # mulcov table:
   mulcov_table = list()
   #
   # avgint table: empty
   avgint_table = list()
   #
   # nslist_table:
   nslist_table = dict()
   # ----------------------------------------------------------------------
   # data table:
   data_table = list()
   #
   # values that are the same for all data points
   row = {
      'integrand':   'Sincidence',
      'hold_out':    False,
      'density':     'gaussian',
      'weight':      '',
      'time_lower':   2000.,
      'time_upper':   2000.,
      'subgroup':     'world',
   }
   for age in [ 0.0 , 100.0 ] :
      for node in [ 'world', 'child_1', 'child_2' ] :
         # Sincidence
         row['meas_value'] = iota_true(age, node)
         row['meas_std']   = row['meas_value'] / 10.0
         row['age_lower']  = age
         row['age_upper']  = age
         row['node']       = node
         data_table.append( copy.copy(row) )
   #
   # ----------------------------------------------------------------------
   # prior_table
   prior_table = [
      { # prior_iota_parent
         'name':     'prior_iota_parent',
         'density':  'uniform',
         'lower':    iota_true(0, 'world') / 10.0,
         'upper':    iota_true(0, 'world') * 10.0,
         'mean':     iota_true(0, 'world') / 2.0 ,
      }, {
         # prior_iota_child
         'name':     'prior_iota_child',
         'density':  'uniform',
         'mean':     0.0 ,
      }
   ]
   # ----------------------------------------------------------------------
   # smooth table
   smooth_table = [
      {  'name':    'smooth_iota_parent',
         'age_id':  [0, 1], # ages 0, 100
         'time_id': [0],
         'fun':     fun_iota_parent
      },{
         'name':    'smooth_iota_child',
         'age_id':  [0],
         'time_id': [0],
         'fun':     fun_iota_child
      }
   ]
   # ----------------------------------------------------------------------
   # rate table:
   rate_table = [
      {  'name':          'iota',
         'parent_smooth': 'smooth_iota_parent',
         'child_smooth':  'smooth_iota_child',
      }
   ]
   # ----------------------------------------------------------------------
   # option_table
   option_table = [
      { 'name':'rate_case',              'value':'iota_pos_rho_zero'   },
      { 'name':'parent_node_name',       'value':'world'               },

      { 'name':'quasi_fixed',            'value':'false'               },
      { 'name':'max_num_iter_fixed',     'value':'50'                  },
      { 'name':'print_level_fixed',      'value':'0'                   },
      { 'name':'tolerance_fixed',        'value':'1e-9'                },

      { 'name':'trace_init_fit_model',   'value':'true'                },
   ]
   # ----------------------------------------------------------------------
   # subgroup_table
   subgroup_table = [ { 'subgroup':'world', 'group':'world' } ]
   # ----------------------------------------------------------------------
   # create database
   dismod_at.create_database(
      file_name,
      age_list,
      time_list,
      integrand_table,
      node_table,
      subgroup_table,
      weight_table,
      covariate_table,
      avgint_table,
      data_table,
      prior_table,
      smooth_table,
      nslist_table,
      rate_table,
      mulcov_table,
      option_table
   )
   # ----------------------------------------------------------------------
   return
# ===========================================================================
# Create database
file_name = 'example.db'
example_db(file_name)
#
program = '../../devel/dismod_at'
dismod_at.system_command_prc([ program, file_name, 'init' ])
#
command = [ program, file_name, 'fit', 'both' ]
stdout  = dismod_at.system_command_prc( command, return_stdout = True )
# -----------------------------------------------------------------------
# Check trace_init_fit_model results
check  = 'Begin dismod_at: fit_model constructor\n'
check += 'Begin cppad_mixed::initialize\n'
check += 'init_ran_like_done_\n'
check += 'init_ran_jac_done_\n'
check += 'init_ran_hes_done_\n'
check += 'init_ldlt_ran_hes_done_\n'
check += 'init_hes_cross_done_\n'
check += 'init_fix_like_done_\n'
check += 'init_fix_con_done_\n'
check += 'End cppad_mixed::initialize\n'
check += 'End dismod_at: fit_model constructor\n'
check += 'Begin cppad_mixed::init_laplace_obj\n'
check += 'init_laplace_obj_fun_done_\n'
check += 'init_laplace_obj_hes_done_\n'
check += 'End cppad_mixed::init_laplace_obj\n'
assert stdout == check
# -----------------------------------------------------------------------
# read database
new                   = False
connection            = dismod_at.create_connection(file_name, new)
node_table            = dismod_at.get_table_dict(connection, 'node')
age_table             = dismod_at.get_table_dict(connection, 'age')
var_table             = dismod_at.get_table_dict(connection, 'var')
fit_var_table         = dismod_at.get_table_dict(connection, 'fit_var')
connection.close()
#
# There are four values for iota in this model
assert len(var_table) == 4
assert len(fit_var_table) == 4
#
# check that the fit is accurate
for var_id in range( len(var_table) ) :
   age_id      = var_table[var_id]['age_id']
   age         = age_table[age_id]['age']
   node_id     = var_table[var_id]['node_id']
   node_name   = node_table[node_id]['node_name']
   true_value  = iota_true(age, node_name)
   if node_name.startswith('child_') :
      parent_value = iota_true(age, 'world')
      true_value   = math.log( true_value / parent_value )
   fit_value   = fit_var_table[var_id]['fit_var_value']
   rel_error   = 1.0 - fit_value/true_value
   assert abs(rel_error) < 1e-6
# -----------------------------------------------------------------------------
print('trace_init.py: OK')
# -----------------------------------------------------------------------------

Input File: example/user/trace_init.py