/*
Programmer: Timothy McNeel, Information Management Services, Inc.
Date: 04/01/2005
This program estimates daily pyramid servings of fruits and vegetables from the 2001
  California Health Interview Survey Adult Data.
Modified 08/29/2007: The variable BASEID was removed because of confidentiality
  concerns.
*/
footnote "Program: /prj/arb/chis/thompson/create.fv.2007_08_29a.sas";

title1 "CHIS 2001 Adult Survey NCI Data File Version 4 (02/15/2005)";

***MODIFY THE FOLLOWING LINES AS APPROPRIATE FOR YOUR LOCAL ENVIRONMENT.;
***AT THE END YOU SHOULD HAVE A DATASET NAMED WORK.CHIS CONTAINING THE 2001 CHIS ADULT RDD DATA.;
  ***Output dataset;
  filename OUTFILE1 pipe 'gzip -c > /prj/arb/chis/thompson/fv.2007_08_29a.v8x.gz';

  ***CHIS 2001 Adult Survey NCI Data File Version 4 (02/15/2005);
  filename CHIS pipe 'gunzip -c /prj/arb/chis/data/2001/ncia4_021505/sasdata/ncia4.v8x.gz';

  ***Formats for CHIS 2001 Adult Survey NCI Data File Version 4 (02/15/2005);
  %include '/prj/arb/chis/data/2001/ncia4_021505/progs/formats.ncia4.sas';

  proc cimport infile=CHIS data=chis;
***MODIFY THE PRECEDING LINES AS APPROPRIATE FOR YOUR LOCAL ENVIRONMENT.;

data chis;
  set chis;

  ***Calculate daily food frequencies and exclude extreme values;
  array aFoodCount[8] ae1-ae8;
  array aFoodUnit[8] ae1unt ae2unt ae3unt ae4unt ae5unt ae6unt ae7unt ae8unt;
  array aFoodPerDay[8] FrtJPD FruitPD FrFryPD OthPotPD DrBeanPD SaladPD OthVegPD TomScPD;
  array aMaxPerDay[8] _temporary_ (4 12 3 3 3 5 9 4);
  do i = 1 to 8;
    if aFoodCount[i]>=0 then do;
      if aFoodCount[i]=0                                      then aFoodPerDay[i]=0;
      else if aFoodUnit[i]=1 and aFoodCount[i]<=aMaxPerDay[i] then aFoodPerDay[i]=aFoodCount[i];
      else if aFoodUnit[i]=2 and aFoodCount[i]<=14            then aFoodPerDay[i]=aFoodCount[i]/7;
      else if aFoodUnit[i]=3 and aFoodCount[i]<=60            then aFoodPerDay[i]=aFoodCount[i]/30;
    end;
  end;

  label AgeGrp = 'Age group';
  if      18 <= srage <= 27 then AgeGrp = 1;
  else if 28 <= srage <= 37 then AgeGrp = 2;
  else if 38 <= srage <= 47 then AgeGrp = 3;
  else if 48 <= srage <= 57 then AgeGrp = 4;
  else if 58 <= srage <= 67 then AgeGrp = 5;
  else if 68 <= srage <= 77 then AgeGrp = 6;
  else if       srage >= 78 then AgeGrp = 7;
run;
/*
proc freq data=chis;
  title2 'Checking new variables';
  tables
    FrtJPD*ae1*ae1unt
    FruitPD*ae2*ae2unt
    FrFryPD*ae3*ae3unt
    OthPotPD*ae4*ae4unt
    DrBeanPD*ae5*ae5unt
    SaladPD*ae6*ae6unt
    OthVegPD*ae7*ae7unt
    TomScPD*ae8*ae8unt
    AgeGrp*srage
    / list missing;
run;
title2;
*/

***Input median portion size in pyramid servings per mention by gender and age for fruits and vegetables analyses, ;
***from the USDA's 1994-96 Continuing Survey of Food Intakes of Individuals (CSFII 94-96) dietary recall data;
data fvps;
  input SRSex AgeGrp FVPSFruit FVPSFrtJ FVPSFrFry FVPSOthPot FVPSSalad FVPSOthVeg FVPSDrBean FVPSTomSc;
  datalines;
1 1 1.301000 2.000000 2.000000 2.000000 0.545000 0.750000 1.374000 0.626000 
1 2 1.301000 1.667500 2.000000 2.000000 0.708000 0.906000 1.047000 0.587000 
1 3 1.229571 1.335000 1.773000 1.999000 0.754500 0.974500 1.065000 0.579000 
1 4 1.227333 1.335000 1.710000 1.999000 0.750000 1.000000 1.227000 0.607000 
1 5 1.168000 1.334000 1.400000 1.914000 0.833500 1.000000 1.000000 0.533600 
1 6 1.168000 1.001000 1.250000 1.544000 0.750000 0.880000 1.000000 0.606667 
1 7 1.052333 1.001000 1.250000 1.508000 0.822500 0.833333 1.114000 0.405000 
2 1 1.168000 1.500500 1.481000 1.544000 0.613500 0.702200 0.964000 0.410000 
2 2 1.168000 1.334000 1.365500 1.544000 0.572500 0.779333 0.684000 0.400000 
2 3 1.168000 1.334000 1.272000 1.528000 0.833333 0.792500 0.800000 0.402000 
2 4 1.168000 1.251250 1.400000 1.544000 1.000000 0.788500 0.687000 0.396000 
2 5 1.150500 1.019500 1.000000 1.499000 0.795500 0.774000 0.822000 0.477500 
2 6 1.083833 1.000500 1.026000 1.516000 0.625000 0.833000 0.807000 0.356000 
2 7 1.000000 1.000500 1.000000 1.272000 0.750000 0.856750 1.000000 0.252000 
;

proc sort data=chis;
  by srsex AgeGrp;
run;

***Add the portion size adjustment variables to the CHIS dataset and estimate daily pyramid servings.;
***Adjust the estimates using regression coefficients from CSFII 94-96.;
data chis;
  merge chis(in=InCHIS) fvps;
  by srsex AgeGrp;
  if InCHIS;

  ***SAS addition operators are used instead of sum statements so that the;
  ***final values will be missing if any of the food variables are missing;

  RawFV = FrtJPD*FVPSFrtJ + FruitPD*FVPSFruit + FrFryPD*FVPSFrFry
    + OthPotPD*FVPSOthPot + DrBeanPD*FVPSDrBean + SaladPD*FVPSSalad
    + OthVegPD*FVPSOthVeg + TomScPD*FVPSTomSc;

  label FV = 'Daily servings of fruits & vegetables (04/01/2005)';
  if      srsex = 1 then FV = ( 0.755105 + 0.815466*(sqrt(RawFV)) )**2;
  else if srsex = 2 then FV = ( 0.711444 + 0.775801*(sqrt(RawFV)) )**2;

  SqRtFV = sqrt(FV);


  RawFVNoFF = FrtJPD*FVPSFrtJ + FruitPD*FVPSFruit
    + OthPotPD*FVPSOthPot + DrBeanPD*FVPSDrBean + SaladPD*FVPSSalad
    + OthVegPD*FVPSOthVeg + TomScPD*FVPSTomSc;

  label FVNoFF = 'Daily servings of fruits & vegetables except French fries (04/01/2005)';
  if      srsex = 1 then FVNoFF = ( 0.771919 + 0.806806*(sqrt(RawFVNoFF)) )**2;
  else if srsex = 2 then FVNoFF = ( 0.694500 + 0.783053*(sqrt(RawFVNoFF)) )**2;

  SqRtFVNoFF = sqrt(FVNoFF);


  RawFVNoFFBns = FrtJPD*FVPSFrtJ + FruitPD*FVPSFruit
    + OthPotPD*FVPSOthPot + SaladPD*FVPSSalad
    + OthVegPD*FVPSOthVeg + TomScPD*FVPSTomSc;

  label FVNoFFBns = 'Daily servings of fruits & vegetables except French fries and beans (04/01/2005)';
  if      srsex = 1 then FVNoFFBns = ( 0.773746 + 0.796528*(sqrt(RawFVNoFFBns)) )**2;
  else if srsex = 2 then FVNoFFBns = ( 0.691176 + 0.776859*(sqrt(RawFVNoFFBns)) )**2;

  SqRtFVNoFFBns = sqrt(FVNoFFBns);


  RawFVNoBns = FrtJPD*FVPSFrtJ + FruitPD*FVPSFruit + FrFryPD*FVPSFrFry
    + OthPotPD*FVPSOthPot + SaladPD*FVPSSalad
    + OthVegPD*FVPSOthVeg + TomScPD*FVPSTomSc;

  label FVNoBns = 'Daily servings of fruits & vegetables except beans (04/01/2005)';
  if      srsex = 1 then FVNoBns = ( 0.729300 + 0.820729*(sqrt(RawFVNoFFBns)) )**2;
  else if srsex = 2 then FVNoBns = ( 0.703858 + 0.772394*(sqrt(RawFVNoFFBns)) )**2;

  SqRtFVNoBns = sqrt(FVNoBns);
run;

proc summary nway data=chis;
  weight rakedw0;
  class srsex;
  var SqRtFV;
  output out=FVmeans mean=SqRtFVMeanBySex;
run;

proc summary nway data=chis;
  weight rakedw0;
  class srsex;
  var SqRtFVNoFF;
  output out=FVNoFFmeans mean=SqRtFVNoFFMeanBySex;
run;

proc summary nway data=chis;
  weight rakedw0;
  class srsex;
  var SqRtFVNoFFBns;
  output out=FVNoFFBnsMeans mean=SqRtFVNoFFBnsMeanBySex;
run;

proc summary nway data=chis;
  weight rakedw0;
  class srsex;
  var SqRtFVNoBns;
  output out=FVNoBnsMeans mean=SqRtFVNoBnsMeanBySex;
run;

proc sort data=chis;           by srsex;
proc sort data=FVmeans;        by srsex;
proc sort data=FVNoFFmeans;    by srsex;
proc sort data=FVNoFFBnsMeans; by srsex;
proc sort data=FVNoBnsMeans;   by srsex;

data chis(keep=puf_id FV FVAdj FVNoFF FVNoFFAdj FVNoFFBns FVNoFFBnsAdj FVNoBns FVNoBnsAdj);
  merge
    chis
    FVmeans(keep=srsex SqRtFVMeanBySex)
    FVNoFFmeans(keep=srsex SqRtFVNoFFMeanBySex)
    FVNoFFBnsMeans(keep=srsex SqRtFVNoFFBnsMeanBySex)
    FVNoBnsMeans(keep=srsex SqRtFVNoBnsMeanBySex)
    ;
  by srsex;

  ***Apply variance adjustment factors from the Eating in America’s Table Study (EATS);

  label FVAdj = 'Variance-adjusted daily servings of fruits & vegetables (04/01/2005)';
  if srsex = 1 then FVAdj = (1.2*(SqRtFV - SqRtFVMeanBySex) + SqRtFVMeanBySex)**2;
  else if srsex = 2 then FVAdj = (1.1*(SqRtFV - SqRtFVMeanBySex) + SqRtFVMeanBySex)**2;

  label FVNoFFAdj = 'Variance-adjusted daily servings of fruits & vegetables except French fries (04/01/2005)';
  if srsex = 1 then FVNoFFAdj = (1.2*(SqRtFVNoFF - SqRtFVNoFFMeanBySex) + SqRtFVNoFFMeanBySex)**2;
  else if srsex = 2 then FVNoFFAdj = (1.1*(SqRtFVNoFF - SqRtFVNoFFMeanBySex) + SqRtFVNoFFMeanBySex)**2;

  label FVNoFFBnsAdj = 'Variance-adjusted daily servings of fruits & vegetables except French fries and beans (04/01/2005)';
  if srsex = 1 then FVNoFFBnsAdj = (1.2*(SqRtFVNoFFBns - SqRtFVNoFFBnsMeanBySex) + SqRtFVNoFFBnsMeanBySex)**2;
  else if srsex = 2 then FVNoFFBnsAdj = (1.1*(SqRtFVNoFFBns - SqRtFVNoFFBnsMeanBySex) + SqRtFVNoFFBnsMeanBySex)**2;

  label FVNoBnsAdj = 'Variance-adjusted daily servings of fruits & vegetables except beans (04/01/2005)';
  if srsex = 1 then FVNoBnsAdj = (1.2*(SqRtFVNoBns - SqRtFVNoBnsMeanBySex) + SqRtFVNoBnsMeanBySex)**2;
  else if srsex = 2 then FVNoBnsAdj = (1.1*(SqRtFVNoBns - SqRtFVNoBnsMeanBySex) + SqRtFVNoBnsMeanBySex)**2;


  array a_vars[*] FV FVAdj FVNoFF FVNoFFAdj FVNoFFBns FVNoFFBnsAdj FVNoBns FVNoBnsAdj;
  do i = 1 to dim(a_vars);
    a_vars[i] = round(a_vars[i],.000001);
  end;
run;

proc sort data=chis;
  by puf_id;
run;

proc contents data=chis;

proc cport data=chis file=OUTFILE1;
