/*****************************************************************************/
/*** The following stored results and files are required for this example: ***/
/***   results from example 1a stored in a SAS IML catalog                 ***/
/***   multivar_distrib_macro_v2.1.sas                                     ***/
/***   hei2010.beanspeas.allocation.macro.sas                              ***/
/***   (see http://riskfactor.cancer.gov/tools/hei/)                       ***/
/***   hei2010.score.macro.v1.2.sas                                        ***/
/***   (see http://riskfactor.cancer.gov/tools/hei/).                      ***/
/***                                                                       ***/
/*** This example includes 5 programs.  Programs 1a, 1b, 1c can be used to ***/
/*** fit a multivariate measurement error model, calculate mean usual      ***/
/*** Healthy Eating Index–2010 (HEI–2010) scores, and perform t-tests      ***/
/*** comparing mean usual HEI–2010 scores using NHANES 2003-2004 data.     ***/
/*** Programs 1d and 1e can be used to extend the analysis to also         ***/
/*** calculate percentiles and cutpoint probablities for both the HEI-2010 ***/
/*** component scores and the HEI-2010 densities, conditional on           ***/
/*** categories of HEI-2010 total scores.                                  ***/
/***                                                                       ***/
/*** This SAS program illustrate the use of the MULTIVAR_DISTRIB macro and ***/
/*** the subsequent calculations of mean usual HEI–2010 scores using       ***/
/*** NHANES 2003-2004 data for smokers, age 20 and older.  This program    ***/
/*** can also be modified to calculate mean usual HEI-2010 scores for      ***/
/*** nonsmokers, age 20 and older.                                         ***/
/***                                                                       ***/
/*** The HEI-2010 nonsmokers and smokers analysis uses balanced repeated   ***/
/*** replication (BRR) variance estimation, so the mean usual HEI-2010     ***/
/*** scores are calculated using the original weights (i.e. weights for    ***/
/*** replicate 0) and using 16 BRR weights (i.e. weights for replicates 1, ***/
/*** 2, ..., 16).  The following replicfirst and repliclast macro          ***/
/*** variables allow specification of a range of replicates.  This example ***/
/*** program considers replicate 0.                                        ***/
/*****************************************************************************/

%let replicfirst = 0;
%let repliclast  = 0;


/***************************************************************/
/*** Specify a title and label for this example and stratum. ***/
/***************************************************************/

title1 "Fit Multivariate Measurement Error Model Using MCMC with 24-Hour Recall as Main Instrument";
title2 "Calculate Mean Usual HEI–2010 Scores and Perform T-Tests Comparing Mean Usual HEI–2010 Scores";
title3 "Calculate Conditional Distributions and Cutpoint Probabilities for HEI-2010 Densities and Scores";
%let example_stratum_label = ex1smok;


/*****************************************************************************/
/*** The path must be modified in the following lines to specify the       ***/
/*** location of the SAS macros and the SAS library that will be used.     ***/
/*** Note that this library is specified using the label provided in the   ***/
/*** above definition of the example_stratum_label macro variable.         ***/
/*****************************************************************************/

%let home = /prj/dcp/statprog/meas.err/develop.public.resources.stat.meth;

*** Include the required macros ***;
%include "&home/include.files.macros/multivar_distrib_macro_v2.1.sas";
%include "&home/include.files.macros/hei2010.beanspeas.allocation.macro.sas";
%include "&home/include.files.macros/hei2010.score.macro.v1.2.sas";

*** Output data library ***;
libname outlib "&home/multivar_surveillance_example1_mcmc_nnc_main24hr/outlib/&example_stratum_label";


proc format ;
  value hei_total_score_gt50_fmt
    0 = '<50'
    1 = '>=50'
    2 = 'All Groups';
run;


/*****************************************************************************/
/*** Use the DISTRIB_HEI_REPLICATE_LOOP macro to call the MULTIVAR_DISTRIB ***/
/*** macro and calculate mean usual HEI-2010 scores for each of the        ***/
/*** replicates specified in the above definition of the replicfirst and   ***/
/*** repliclast macro variables.                                           ***/
/*****************************************************************************/


/*****************************************************************************/
/*                                                                           */
/* Description of the MULTIVAR_DISTRIB macro                                 */
/*                                                                           */
/*****************************************************************************/
/*                                                                           */
/*                                                                           */
/* The MULTIVAR_DISTRIB macro uses parameter estimates from the              */
/* MULTIVAR_MCMC macro to generate a multivariate Monte Carlo distribution   */
/* of the usual intakes for the dietary components specified in the          */
/* multivariate measurement error model fit using the MULTIVAR_MCMC macro.   */
/* The MULTIVAR_DISTRIB macro can also use values that were stored using     */
/* the macro parameter "optional_iml_store_names" of the MULTIVAR_MCMC       */
/* macro.  The MULTIVAR_DISTRIB macro allows specification of 1 or 2 lists   */
/* of covariates.  If 2 lists of covariates are specified, the usual intake  */
/* is calculated as the weighted average of the usual intake calculated      */
/* using covariate list 1 and the usual intake calculated using covariate    */
/* list 2.  This weighted average is calculated using the values specified   */
/* for the "set_value_for_weight_cov_list1" and                              */
/* "set_value_for_weight_cov_list2" macro parameters.                        */
/*                                                                           */
/* References:                                                               */
/*                                                                           */
/*   Zhang S, Krebs-Smith SM, Midthune D, Perez A, Buckman DW, Kipnis V,     */
/*   Freedman LS, Dodd KW, Carroll RJ. Fitting a bivariate measurement error */
/*   model for episodically consumed dietary components. Int J Biostat       */
/*   2011;7(1):Article 1.                                                    */
/*                                                                           */
/*   Zhang S, Midthune D, Guenther PM, Krebs-Smith SM, Kipnis V, Dodd KW,    */
/*   Buckman DW, Tooze JA, Freedman L, Carroll RJ. A new multivariate        */
/*   measurement error model with zero-inflated dietary data, and its        */
/*   application to dietary assessment. Ann Appl Stat 2011 Jun;5(2B):        */
/*   1456-87.                                                                */
/*                                                                           */
/*                                                                           */
/* The syntax for calling the MULTIVAR_DISTRIB macro is:                     */
/*                                                                           */
/* %multivar_distrib(multivar_mcmc_out_lib           =,                      */
/*                   multivar_mcmc_out_store_label   =,                      */
/*                   t_weightavg_covariates_list1    =,                      */
/*                   t_weightavg_covariates_list2    =,                      */
/*                   set_value_for_weight_cov_list1  =,                      */
/*                   set_value_for_weight_cov_list2  =,                      */
/*                   nev_consumers_epis1             =,                      */
/*                   covars_prob_consumer_epis1      =,                      */
/*                   optional_input_data             =,                      */
/*                   optional_input_data_var_list    =,                      */
/*                   optional_input_mcmc_u_out_data  =,                      */
/*                   additional_output_var_list      =,                      */
/*                   additional_output_subject_var   =,                      */
/*                   output_mcmc_weight_var          =,                      */
/*                   set_seed_distrib                =,                      */
/*                   set_number_monte_carlo_rand_obs =,                      */
/*                   print                           =                       */
/*                   );                                                      */
/*                                                                           */
/*  where                                                                    */
/*                                                                           */
/*  "multivar_mcmc_out_lib"            Specifies the same SAS library that   */
/*                                     was specified for the "out_lib" macro */
/*                                     parameter of the MULTIVAR_MCMC macro. */
/*                                     The MULTIVAR_DISTRIB macro will use   */
/*                                     the parameter estimates and can use   */
/*                                     covariates that were stored by the    */
/*                                     MULTIVAR_MCMC macro.  Covariates can  */
/*                                     be stored using the                   */
/*                                     "optional_iml_store_names" macro      */
/*                                     parameter of the MULTIVAR_MCMC macro. */
/*                                                                           */
/*  "multivar_mcmc_out_store_label"    Specifies the same label that was     */
/*                                     specified for the "out_store_label"   */
/*                                     macro parameter of the MULTIVAR_MCMC  */
/*                                     macro.  The MULTIVAR_DISTRIB macro    */
/*                                     will use the parameter estimates and  */
/*                                     can use covariates that were stored   */
/*                                     by the MULTIVAR_MCMC macro.           */
/*                                     Covariates can be stored using the    */
/*                                     "optional_iml_store_names" macro      */
/*                                     parameter of the MULTIVAR_MCMC macro. */
/*                                                                           */
/*  "t_weightavg_covariates_list1"     Specifies a list of covariates that   */
/*                                     correspond to the beta parameter      */
/*                                     estimates from the multivariate       */
/*                                     measurement error model fit using the */
/*                                     MULTIVAR_MCMC macro.  The             */
/*                                     calculations in the MULTIVAR_DISTRIB  */
/*                                     macro are based on the assumption     */
/*                                     that the multivariate measurement     */
/*                                     error model was fit using the same    */
/*                                     set of covariates for each part of    */
/*                                     the model.  The MULTIVAR_DISTRIB      */
/*                                     macro creates "constant0" which can   */
/*                                     be used as a variable that equals 0   */
/*                                     for each observation, and the macro   */
/*                                     creates "constant1" which can be used */
/*                                     as a variable that equals 1 for each  */
/*                                     observation.  "Constant0" and         */
/*                                     "constant1" can be used when          */
/*                                     specifying covariates of interest.    */
/*                                     If the optional macro parameter       */
/*                                     "t_weightavg_covariates_list2" is     */
/*                                     used to specify a second list of      */
/*                                     covariates, then usual intake is      */
/*                                     calculated as the weighted average of */
/*                                     the usual intake calculated using     */
/*                                     covariate list 1 and the usual intake */
/*                                     calculated using covariate list 2.    */
/*                                     This weighted average is calculated   */
/*                                     using the values specified for the    */
/*                                     "set_value_for_weight_cov_list1" and  */
/*                                     "set_value_for_weight_cov_list2"      */
/*                                     macro parameters.  If the second list */
/*                                     of covariates is not specified, then  */
/*                                     the Monte Carlo distribution of usual */
/*                                     intakes is only based on the          */
/*                                     covariates specified for the          */
/*                                     "t_weightavg_covariates_list1" macro  */
/*                                     parameter.                            */
/*                                                                           */
/*  "t_weightavg_covariates_list2"     Specifies an optional second list of  */
/*                                     covariates.  If this list of          */
/*                                     covariates is specified, the usual    */
/*                                     intake is calculated as the weighted  */
/*                                     average of the usual intake           */
/*                                     calculated using covariate list 1 and */
/*                                     the usual intake calculated using     */
/*                                     covariate list 2.  This weighted      */
/*                                     average is calculated using the       */
/*                                     values specified for the              */
/*                                     "set_value_for_weight_cov_list1" and  */
/*                                     "set_value_for_weight_cov_list2"      */
/*                                     macro parameters.                     */
/*                                                                           */
/*  "set_value_for_weight_cov_list1"   Specifies a value that is used to     */
/*                                     calculate a weighted average when 2   */
/*                                     covariate lists are provided as       */
/*                                     described for the                     */
/*                                     "t_weightavg_covariates_list1" and    */
/*                                     "t_weightavg_covariates_list2" macro  */
/*                                     parameters.  The value specified is   */
/*                                     used as the numerator of the weight.  */
/*                                     The denominator of the weight is      */
/*                                     defined as the sum of the values from */
/*                                     the "set_value_for_weight_cov_list1"  */
/*                                     and "set_value_for_weight_cov_list2"  */
/*                                     macro parameters.  The default value  */
/*                                     is 1 which corresponds to the         */
/*                                     situation where the                   */
/*                                     "t_weightavg_covariates_list1" macro  */
/*                                     parameter is used to provide          */
/*                                     covariate list 1 and the              */
/*                                     "t_weightavg_covariates_list2" macro  */
/*                                     parameter is not used.                */
/*                                                                           */
/*  "set_value_for_weight_cov_list2"   Specifies a value that is used to     */
/*                                     calculate a weighted average when 2   */
/*                                     covariate lists are provided as       */
/*                                     described for the                     */
/*                                     "t_weightavg_covariates_list1" and    */
/*                                     "t_weightavg_covariates_list2" macro  */
/*                                     parameters.  The value specified is   */
/*                                     used as the numerator of the weight.  */
/*                                     The denominator of the weight is      */
/*                                     defined as the sum of the values from */
/*                                     the "set_value_for_weight_cov_list1"  */
/*                                     and "set_value_for_weight_cov_list2"  */
/*                                     macro parameters.  The default value  */
/*                                     is 0 which corresponds to the         */
/*                                     situation where the                   */
/*                                     "t_weightavg_covariates_list2" macro  */
/*                                     parameter is not used.                */
/*                                                                           */
/*  "nev_consumers_epis1"              Specifies the same setting that was   */
/*                                     specified for the                     */
/*                                     "nev_consumers_epis1" macro parameter */
/*                                     of the MULTIVAR_MCMC macro.  The      */
/*                                     default value is "n".                 */
/*                                                                           */
/*  "covars_prob_consumer_epis1"       Specifies the same list of covariates */
/*                                     that was specified for the            */
/*                                     "covars_prob_consumer_epis1"          */
/*                                     macro parameter of the MULTIVAR_MCMC  */
/*                                     macro.                                */
/*                                                                           */
/*  "optional_input_data"              Specifies an optional SAS data set    */
/*                                     that provides the additional input    */
/*                                     variables listed in the               */
/*                                     "optional_input_data_var_list" macro  */
/*                                     parameter.  If this data set is       */
/*                                     specified and the                     */
/*                                     "additional_output_subject_var" macro */
/*                                     parameter specifies a subject         */
/*                                     identification variable, then the     */
/*                                     subject identification variable is    */
/*                                     included in the output Monte Carlo    */
/*                                     data set.  When an optional input     */
/*                                     data set is not specified, the Monte  */
/*                                     Carlo "mc_t_distrib_out" data set     */
/*                                     includes the "weight_nw_sumw"         */
/*                                     variable created and stored by the    */
/*                                     MULTIVAR_MCMC macro (see the          */
/*                                     MULTIVAR_MCMC "weight_var" macro      */
/*                                     parameter for details).  When an      */
/*                                     optional input data set is specified, */
/*                                     an advanced user can use the          */
/*                                     "output_mcmc_weight_var" macro        */
/*                                     parameter to include the              */
/*                                     "weight_nw_sumw" variable in the      */
/*                                     "mc_t_distrib_out" data set if        */
/*                                     appropriate.  When an optional input  */
/*                                     data set and both the                 */
/*                                     "additional_output_var_list" and      */
/*                                     "additional_output_subject_var" macro */
/*                                     parameters are specified, an advanced */
/*                                     user can include these additional     */
/*                                     variables in the "mc_t_distrib_out"   */
/*                                     data set.                             */
/*                                                                           */
/*  "optional_input_data_var_list"     Specifies optional input variables.   */
/*                                     For details, see the description for  */
/*                                     the "optional_input_data" macro       */
/*                                     parameter.  This variable list should */
/*                                     not include names of variables stored */
/*                                     by the MULTIVAR_MCMC macro since the  */
/*                                     stored values are brought into the    */
/*                                     macro as described for the            */
/*                                     "multivar_mcmc_out_lib" and           */
/*                                     "multivar_mcmc_out_store_label" macro */
/*                                     parameters.                           */
/*                                                                           */
/*  "optional_input_mcmc_u_out_data"   Specifies an optional SAS data set    */
/*                                     that provides u matrices saved from   */
/*                                     iterations of the MULTIVAR_MCMC       */
/*                                     macro.  Considering the MULTIVAR_MCMC */
/*                                     "set_number_post_mcmc_u_out" macro    */
/*                                     parameter, the variables required for */
/*                                     this input data set include:          */
/*                                     - "Iteration" which indexes the first */
/*                                       u matrix through last u matrix in   */
/*                                       the data set,                       */
/*                                     - the variable specified by the       */
/*                                       "subject" macro parameter from the  */
/*                                       MULTIVAR_MCMC macro,                */
/*                                     - u_col1, u_col2, ..., u_colL where L */
/*                                       is the number of columns of the u   */
/*                                       matrix.  L is defined in the        */
/*                                       description of the                  */
/*                                       "out_store_label" macro parameter   */
/*                                       of the MULTIVAR_MCMC macro.         */
/*                                     The MULTIVAR_DISTRIB macro will use   */
/*                                     the first (i.e. iteration 1) u matrix */
/*                                     under the assumption that the rows of */
/*                                     this matrix are still sorted in the   */
/*                                     same order as when the MULTIVAR_MCMC  */
/*                                     macro generated this u matrix, and    */
/*                                     the MULTIVAR_DISTRIB macro will       */
/*                                     generate and append additional rows   */
/*                                     to this matrix as needed to obtain    */
/*                                     the dimensions as specified and       */
/*                                     implied by other user input to this   */
/*                                     macro.  The MULTIVAR_DISTRIB macro    */
/*                                     will use this same approach for each  */
/*                                     subsequent u matrix and iteration.    */
/*                                                                           */
/*  "additional_output_var_list"       Specifies additional output           */
/*                                     variables.  For details, see the      */
/*                                     description for the                   */
/*                                     "optional_input_data" macro           */
/*                                     parameter.                            */
/*                                                                           */
/*  "additional_output_subject_var"    Specifies a subject identification    */
/*                                     variable to be included in the output */
/*                                     Monte Carlo data set.  For details,   */
/*                                     see the description for the           */
/*                                     "optional_input_data" macro           */
/*                                     parameter.                            */
/*                                                                           */
/*  "output_mcmc_weight_var"           When an optional input data set is    */
/*                                     specified, an advanced user can       */
/*                                     specify "output_mcmc_weight_var=y" or */
/*                                     "output_mcmc_weight_var=Y" to include */
/*                                     the "weight_nw_sumw" variable in the  */
/*                                     "mc_t_distrib_out" data set if        */
/*                                     appropriate.  For details, see the    */
/*                                     description for the                   */
/*                                     "optional_input_data" macro           */
/*                                     parameter.  The default value is "n". */
/*                                                                           */
/*  "set_seed_distrib"                 Specifies a seed for random number    */
/*                                     generation for the Monte Carlo method */
/*                                     used in the macro.  If                */
/*                                     "set_seed_distrib" is not specified,  */
/*                                     the SAS IML procedure will generate   */
/*                                     an initial seed value from the system */
/*                                     clock.                                */
/*                                                                           */
/*  "set_number_monte_carlo_rand_obs"  Specifies the number of observations  */
/*                                     of the random effects vector to       */
/*                                     generate for an individual.  If an    */
/*                                     "optional_input_mcmc_u_out_data" SAS  */
/*                                     data set is specified, then the       */
/*                                     "set_number_monte_carlo_rand_obs"     */
/*                                     macro parameter is not used, and the  */
/*                                     number of observations of the random  */
/*                                     effects vector is determined by the   */
/*                                     maximum value of the "Iteration"      */
/*                                     index variable.  Additional details   */
/*                                     are provided in the description of    */
/*                                     the "optional_input_mcmc_u_out_data"  */
/*                                     macro parameter.  If an               */
/*                                     "optional_input_mcmc_u_out_data" SAS  */
/*                                     data set is not specified, the        */
/*                                     number of observations in the Monte   */
/*                                     Carlo data set is the product of the  */
/*                                     "number of individuals" and the       */
/*                                     specified value for the               */
/*                                     "set_number_monte_carlo_rand_obs"     */
/*                                     macro parameter.  The default value   */
/*                                     is 500.                               */
/*                                                                           */
/*  "print"                            If "print=n" or "print=N" then macro  */
/*                                     information is not printed.           */
/*                                     Otherwise the macro prints            */
/*                                     information regarding the IML storage */
/*                                     catalog used and prints some stored   */
/*                                     and specified values.  The default    */
/*                                     value is "y".                         */
/*                                                                           */
/*****************************************************************************/
/*                                                                           */
/* Macro Output Data Set:  mc_t_distrib_out                                  */
/*                                                                           */
/* The MULTIVAR_DISTRIB macro produces a SAS data set "mc_t_distrib_out"     */
/* that includes usual intake variables named:                               */
/*          mc_t1   mc_t2   ...   mc_tp                                      */
/* where p = num_epis_diet_comp + num_daily_diet_comp, and the dietary       */
/* component order (1, 2, ..., p) is equivalent to the order used when the   */
/* multivariate measurement error model was fit using the MULTIVAR_MCMC      */
/* macro.  When the "t_weightavg_covariates_list2" macro parameter is not    */
/* used, the output data set "mc_t_distrib_out" includes the following       */
/* variables used to calculate the usual intake variables:                   */
/*          mc_prob1   mc_prob2   ...   mc_probq                             */
/*          mc_amount1   mc_amount2   ...   mc_amountq                       */
/* where q = num_epis_diet_comp, and these variables correspond to the       */
/* probability to consume and the intake amount on consumption days.         */
/*                                                                           */
/* When the "optional_input_data" macro parameter is not used to specify an  */
/* optional input data set, the "mc_t_distrib_out" data set also             */
/* includes the variable:                                                    */
/*          weight_nw_sumw                                                   */
/* created and stored by the MULTIVAR_MCMC macro (see the MULTIVAR_MCMC      */
/* "weight_var" macro parameter for details).                                */
/*                                                                           */
/* Additional variables can be stored in the "mc_t_distrib_out" data set     */
/* (see the macro parameter descriptions for details).                       */
/*                                                                           */
/* If an "optional_input_mcmc_u_out_data" SAS data set is not specified, the */
/* number of observations in the output data set is the product of the       */
/* "number of individuals" and the specified value for the                   */
/* "set_number_monte_carlo_rand_obs" macro parameter (see the macro          */
/* parameter descriptions for details).                                      */
/*                                                                           */
/*****************************************************************************/
/*****************************************************************************/
/*****************************************************************************/



%macro distrib_hei_replicate_loop;

  %do replicnum = &replicfirst %to &repliclast;


    %let outlabel = &example_stratum_label._replic&replicnum;
    %let loopseed_distrib = %eval(89009890 + &replicnum * 10000);


    /*************************************************************************/
    /*** Call the MULTIVAR_DISTRIB macro where covariate list 1 and        ***/
    /*** covariate list 2 are specified for week days and weekends,        ***/
    /*** respectively.  Note that the multivariate measurement error model ***/
    /*** was fit using the MULTIVAR_MCMC macro with the following          ***/
    /*** covariates:                                                       ***/
    /***   constant1, recallday2, weekend, and other covariates            ***/
    /***   (see example 1a).                                               ***/
    /*************************************************************************/;

    %multivar_distrib(multivar_mcmc_out_lib           = outlib,
                      multivar_mcmc_out_store_label   = surv_mcmc_nc_main24hr&replicnum,
                      t_weightavg_covariates_list1    = constant1 constant0 constant0 age20_30 age31_50 age51_70 raceeth_mexam
                                                        raceeth_nhsblk male male_age20_30 male_age31_50 male_age51_70
                                                        male_raceeth_mexam male_raceeth_nhsblk,
                      t_weightavg_covariates_list2    = constant1 constant0 constant1 age20_30 age31_50 age51_70 raceeth_mexam
                                                        raceeth_nhsblk male male_age20_30 male_age31_50 male_age51_70
                                                        male_raceeth_mexam male_raceeth_nhsblk,
                      set_value_for_weight_cov_list1  = 4,
                      set_value_for_weight_cov_list2  = 3,
                      nev_consumers_epis1             = ,
                      covars_prob_consumer_epis1      = ,
                      optional_input_data             = ,
                      optional_input_data_var_list    = ,
                      optional_input_mcmc_u_out_data  = ,
                      additional_output_var_list      = ,
                      additional_output_subject_var   = replicaterowid,
                      output_mcmc_weight_var          = ,
                      set_seed_distrib                = &loopseed_distrib,
                      set_number_monte_carlo_rand_obs = 500,
                      print                           =
                      );


    data usintake_mc_t;
      set mc_t_distrib_out;


      /***********************************************************************/
      /*** Divide weights by set_number_monte_carlo_rand_obs value used in ***/
      /*** the MULTIVAR_DISTRIB macro.                                     ***/
      /***********************************************************************/

      weight_nw_sumw_div500 = weight_nw_sumw / 500;


      /********************************************************/
      /*** Assign variable names for usual intakes from the ***/
      /*** MULTIVAR_DISTRIB macro.                          ***/
      /********************************************************/

      usintake_frtjuice        = mc_t1;
      usintake_wholefrt        = mc_t2;
      usintake_obeangrn        = mc_t3;
      usintake_olegumes        = mc_t4;
      usintake_g_whl           = mc_t5;
      usintake_oseaplant       = mc_t6;
      usintake_ononbeangrn     = mc_t7;
      usintake_g_nwhl          = mc_t8;
      usintake_d_total         = mc_t9;
      usintake_ononseaplant    = mc_t10;
      usintake_tsodi           = mc_t11;
      usintake_tmfat           = mc_t12;
      usintake_tpfat           = mc_t13;
      usintake_tsfat           = mc_t14;
      usintake_emptycal        = mc_t15;
      usintake_nonemptycal     = mc_t16;


      usintake_tot_fruit       = usintake_frtjuice + usintake_wholefrt;
      usintake_tot_energy      = usintake_emptycal + usintake_nonemptycal;


      usintake_otot_veg        = usintake_obeangrn + usintake_ononbeangrn;
      usintake_otot_prot       = usintake_oseaplant + usintake_ononseaplant;


      output usintake_mc_t;
    run;


    /***************************************************************************/
    /*** Allocate beans and peas (i.e. legumes) using a SAS macro, LEG2010A. ***/
    /***************************************************************************/

    %leg2010a(indat    = usintake_mc_t,
              kcal     = usintake_tot_energy,
              allmeat  = usintake_otot_prot,
              seaplant = usintake_oseaplant,
              v_total  = usintake_otot_veg,
              v_drkgr  = usintake_obeangrn,
              legumes  = usintake_olegumes,
              outdat   = usintake_leg_alloc
              );


    data usintake_leg_alloc;
      set usintake_leg_alloc;

      rename legume_added_usintake_otot_prot = usintake_tot_prot
             legume_added_usintake_oseaplant = usintake_seaplant
             legume_added_usintake_otot_veg  = usintake_tot_veg
             legume_added_beangrn            = usintake_beangrn
             ;

      usintake_monopoly = usintake_tmfat + usintake_tpfat;
    run;


    /*********************************************************************/
    /*** Calculate usual HEI-2010 scores using the multivariate Monte  ***/
    /*** Carlo distribution of usual intakes and a SAS macro, HEI2010. ***/
    /*********************************************************************/

    %hei2010(indat      = usintake_leg_alloc,
             kcal       = usintake_tot_energy,
             f_total    = usintake_tot_fruit,
             wholefrt   = usintake_wholefrt,
             lv_total   = usintake_tot_veg,
             lbeangrn   = usintake_beangrn,
             g_whl      = usintake_g_whl,
             d_total    = usintake_d_total,
             lallmeat   = usintake_tot_prot,
             lseaplant  = usintake_seaplant,
             monopoly   = usintake_monopoly,
             sfat       = usintake_tsfat,
             g_nwhl     = usintake_g_nwhl,
             sodi       = usintake_tsodi,
             emptycal10 = usintake_emptycal,
             outdat     = usintake_dens_hei
             );


    /*********************************************/
    /*** Calculate mean usual HEI-2010 scores. ***/
    /*********************************************/

    proc means data=usintake_dens_hei mean noprint;
      weight weight_nw_sumw_div500;
      var heix1_totalveg
          heix2_green_and_bean
          heix3_totalfruit
          heix4_wholefruit
          heix5_wholegrain
          heix6_totaldairy
          heix7_totprot
          heix8_seaplant_prot
          heix9_fattyacid
          heix10_sodium
          heix11_refinedgrain
          heix12_sofaas
          hei2010_total_score
          ;
      output out  = hei_mean_out1rec
             mean = mean_smok_hei_tot_veg
                    mean_smok_hei_beangrn
                    mean_smok_hei_tot_fruit
                    mean_smok_hei_wholefrt
                    mean_smok_hei_g_whl
                    mean_smok_hei_d_total
                    mean_smok_hei_tot_prot
                    mean_smok_hei_seaplant
                    mean_smok_hei_fat_acids_ratio
                    mean_smok_hei_tsodi
                    mean_smok_hei_g_nwhl
                    mean_smok_hei_emptycal
                    mean_smok_hei_total
                    ;
    run;


    proc print data=hei_mean_out1rec;
      var mean_smok_hei_tot_fruit
          mean_smok_hei_wholefrt
          mean_smok_hei_tot_veg
          mean_smok_hei_beangrn
          mean_smok_hei_g_whl
          mean_smok_hei_d_total
          mean_smok_hei_tot_prot
          mean_smok_hei_seaplant
          mean_smok_hei_fat_acids_ratio
          mean_smok_hei_g_nwhl
          mean_smok_hei_tsodi
          mean_smok_hei_emptycal
          mean_smok_hei_total
          ;
    run;


  %end; *** end loop for replicates ***;


%mend;



%distrib_hei_replicate_loop;
