\(\newcommand{\B}[1]{ {\bf #1} }\) \(\newcommand{\R}[1]{ {\rm #1} }\)
continue_cascade.py¶
View page sourcecontinue_cascade: Python Source Code¶
# ----------------------------------------------------------------------------
# imports
# ----------------------------------------------------------------------------
import math
import sys
import os
import copy
import time
import csv
import shutil
import dismod_at
from math import exp
#
# import at_cascade with a preference current directory version
current_directory = os.getcwd()
if os.path.isfile( current_directory + '/at_cascade/__init__.py' ) :
sys.path.insert(0, current_directory)
import at_cascade
# -----------------------------------------------------------------------------
# global variables
# -----------------------------------------------------------------------------
# BEGIN fit_goal_set
first_fit_goal_set = { 'n3', 'n4', 'n2' }
second_fit_goal_set = { 'n5', 'n6' }
# END fit_goal_set
# BEGIN option_all
option_all = {
'result_dir': 'build/example',
'root_node_name': 'n0',
'max_number_cpu': '2',
}
option_all['root_database'] = option_all['result_dir'] + '/root.db'
# END option_all
# ----------------------------------------------------------------------------
# functions
# ----------------------------------------------------------------------------
# BEGIN rate_true
def rate_true(rate, a, t, n, c) :
iota_true = {
'n3' : 0.04,
'n4' : 0.05,
'n5' : 0.06,
'n6' : 0.07,
}
iota_true['n1'] = (iota_true['n3'] + iota_true['n4']) / 2.0
iota_true['n2'] = (iota_true['n5'] + iota_true['n6']) / 2.0
iota_true['n0'] = (iota_true['n1'] + iota_true['n2']) / 2.0
if rate == 'iota' :
return iota_true[n]
return 0.0
# END rate_true
# ----------------------------------------------------------------------------
def root_node_db(file_name) :
#
# BEGIN iota_mean
iota_mean = rate_true('iota', None, None, 'n0', None)
# END iota_mean
#
# prior_table
prior_table = list()
prior_table.append(
# BEGIN parent_value_prior
{ 'name': 'parent_value_prior',
'density': 'gaussian',
'lower': iota_mean / 10.0,
'upper': iota_mean * 10.0,
'mean': iota_mean,
'std': iota_mean,
'eta': iota_mean * 1e-3
}
# END parent_value_prior
)
prior_table.append(
# BEGIN child_value_prior
{ 'name': 'child_value_prior',
'density': 'gaussian',
'mean': 0.0,
'std': 10.0,
}
# END child_value_prior
)
#
# smooth_table
smooth_table = list()
#
# parent_smooth
fun = lambda a, t : ('parent_value_prior', None, None)
smooth_table.append({
'name': 'parent_smooth',
'age_id': [0],
'time_id': [0],
'fun': fun,
})
#
# child_smooth
fun = lambda a, t : ('child_value_prior', None, None)
smooth_table.append({
'name': 'child_smooth',
'age_id': [0],
'time_id': [0],
'fun': fun,
})
#
# node_table
node_table = [
{ 'name':'n0', 'parent':'' },
{ 'name':'n1', 'parent':'n0' },
{ 'name':'n2', 'parent':'n0' },
{ 'name':'n3', 'parent':'n1' },
{ 'name':'n4', 'parent':'n1' },
{ 'name':'n5', 'parent':'n2' },
{ 'name':'n6', 'parent':'n2' },
]
#
# rate_table
rate_table = [ {
'name': 'iota',
'parent_smooth': 'parent_smooth',
'child_smooth': 'child_smooth',
} ]
#
# covariate_table
covariate_table = list()
#
# mulcov_table
mulcov_table = list()
#
# subgroup_table
subgroup_table = [ {'subgroup': 'world', 'group':'world'} ]
#
# integrand_table
integrand_table = [ {'name':'Sincidence'} ]
#
# avgint_table
avgint_table = list()
row = {
'node': 'n0',
'subgroup': 'world',
'weight': '',
'age_lower': 50.0,
'age_upper': 50.0,
'time_lower': 2000.0,
'time_upper': 2000.0,
'integrand': 'Sincidence',
}
avgint_table.append( copy.copy(row) )
#
# data_table
data_table = list()
leaf_set = { 'n3', 'n4', 'n5', 'n6' }
row = {
'subgroup': 'world',
'weight': '',
'age_lower': 50.0,
'age_upper': 50.0,
'time_lower': 2000.0,
'time_upper': 2000.0,
'integrand': 'Sincidence',
'density': 'gaussian',
'hold_out': False,
}
for node in leaf_set :
meas_value = rate_true('iota', None, None, node, None)
row['node'] = node
row['meas_value'] = meas_value
row['meas_std'] = meas_value / 10.0
data_table.append( copy.copy(row) )
#
# age_grid
age_grid = [ 0.0, 100.0 ]
#
# time_grid
time_grid = [ 2000.0 ]
#
# weight table:
weight_table = list()
#
# nslist_table
nslist_table = dict()
#
# option_table
# print_level_fixed is 5 and max_number_cpu > 1 so optimizer trace
# is printed to a file in same directory as corresponding database.
option_table = [
{ 'name':'parent_node_name', 'value':'n0'},
{ 'name':'rate_case', 'value':'iota_pos_rho_zero'},
{ 'name': 'zero_sum_child_rate', 'value':'iota'},
{ 'name':'quasi_fixed', 'value':'false'},
{ 'name':'max_num_iter_fixed', 'value':'50'},
{ 'name':'tolerance_fixed', 'value':'1e-8'},
{ 'name':'print_level_fixed', 'value':'5'},
]
# ----------------------------------------------------------------------
# create database
dismod_at.create_database(
file_name,
age_grid,
time_grid,
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
)
# ----------------------------------------------------------------------------
def two_fit_goal_set_example(result_dir, avgint_table) :
# -------------------------------------------------------------------------
#
# root_database
root_database = option_all['root_database']
#
# all_node_database
all_node_database = f'{result_dir}/all_node.db'
#
# root_fit_dir
root_fit_dir = f'{result_dir}/n0'
at_cascade.empty_directory(root_fit_dir)
#
# cascade starting at n0
at_cascade.cascade_root_node(
all_node_database = all_node_database ,
fit_goal_set = first_fit_goal_set ,
)
#
# continue starting at at n2
fit_database = root_fit_dir + '/n2/dismod.db'
at_cascade.continue_cascade(
all_node_database = all_node_database ,
fit_database = fit_database ,
fit_goal_set = second_fit_goal_set ,
)
#
# check leaf node results
leaf_dir_list = [ 'n0/n1/n3', 'n0/n1/n4', 'n0/n2/n5', 'n0/n2/n6' ]
for leaf_dir in leaf_dir_list :
leaf_database = f'{result_dir}/{leaf_dir}/dismod.db'
at_cascade.check_cascade_node(
rate_true = rate_true,
all_node_database = all_node_database,
fit_database = leaf_database,
avgint_table = avgint_table,
relative_tolerance = 1e-7,
)
# ----------------------------------------------------------------------------
def one_fit_goal_set_example(result_dir, avgint_table) :
#
# root_database
root_database = option_all['root_database']
#
# all_node_database
all_node_database = f'{result_dir}/all_node.db'
#
# root_fit_dir
root_fit_dir = f'{result_dir}/n0'
at_cascade.empty_directory(root_fit_dir)
#
# cascade starting at n0
fit_goal_set = first_fit_goal_set | second_fit_goal_set
at_cascade.cascade_root_node(
all_node_database = all_node_database ,
fit_goal_set = fit_goal_set ,
)
#
# check leaf node results
leaf_dir_list = [ 'n0/n1/n3', 'n0/n1/n4', 'n0/n2/n5', 'n0/n2/n6' ]
for leaf_dir in leaf_dir_list :
leaf_database = f'{result_dir}/{leaf_dir}/dismod.db'
at_cascade.check_cascade_node(
rate_true = rate_true,
all_node_database = all_node_database,
fit_database = leaf_database,
avgint_table = avgint_table,
relative_tolerance = 1e-7,
)
# ----------------------------------------------------------------------------
# main
# ----------------------------------------------------------------------------
def main() :
# -------------------------------------------------------------------------
# result_dir
result_dir = option_all['result_dir']
at_cascade.empty_directory(result_dir)
#
# root_database
root_database = option_all['root_database']
root_node_db(root_database)
#
# avgint_table
# This also erases the avgint table from root_database
avgint_table = at_cascade.extract_avgint( root_database )
#
# all_node_database
all_node_database = f'{result_dir}/all_node.db'
at_cascade.create_all_node_db(
all_node_database = all_node_database,
option_all = option_all,
)
#
# example using continue_cascade
two_fit_goal_set_example(result_dir, avgint_table)
#
# same calculation without continue_cascade
one_fit_goal_set_example(result_dir, avgint_table)
#
# Without this, the mac will try to execute main on each processor.
if __name__ == '__main__' :
main()
print('max_fit_option: OK')