Manual: simpen

11/18/05

1.  Introduction
2.  Examples
    a.  Generate model only
    b.  Generate model and data
    c.  Generate data only
3.  Input files
4.  Model generating fitness function 
5.  Genotype data simulation
6.  Screen output 
7.  File output


1.  Introduction

The program simpen combines a penetrance model generator with a genotype
data simulator.  Simpen can be used to evolve a new penetrance model only or
it can use that penetrance model to generate the simulated genotype data as
specified in the configuration file.  Alternatively, simpen can accept a 
previously defined model and generate the genotype data.

The program utilizes a genetic algorithm to evolve a penetrance
model meeting the specifications of the user.  The program minimizes the 
marginal penetrance variance to simulate a model with minimal main effects
while also optimizing heritability, table variance, and average marginal penetrance
as selected by the user.

Simpen is linked to a data simulator that will use the model evolved by simpen
to create data files of individuals with genotypes.  This program also accepts 
a datasim file that lists the parameters for running the simulator.  
The examples provide additional information on these options.


2. Examples

    Simpen takes 2 arguments.  The first is the name of the configuration file
    and the second is the base name to use in creating output files.  Sample 
    files are at the end of this document.  Simpen automatically determines
    whether the input file is a file for producing models or whether the file
    specifies the production of data only.

a.  Running simpen to only produce a model meeting specified requirements:

    Command:
    simpen modelonly.simpen modelonly
    
    Input file parameters:

    SEED        2000
    LOCI        2
    FREQ        0.8 0.2
    TABLEVAR    0.10
    USETABLEVAR OFF
    HERIT       0.03
    HERITWEIGHT 2
    ODDSRATIO   5
    ODDSWEIGHT  .02
    MARGVAR     0.00001
    MARGWEIGHT  8
    MODELS      1
    GEN         3000
    POPSIZE     100
    DEMES       10
    MUTATE      0.01
    CROSS       0.6
    VERBOSE     ON
    SUBMODELS   ON
    UPDATE      100
    SIMTYPE     prob
    
    
    Output files:
    modelonly.out --> Lists the penetrance of the models evolved and the stats.
    modelonly.1.datasim -->  File used with data simulator.  The '1' refers to
                             this using the first model evolved by simpen.
    modelonly.1.smod    -->  File used with data simulator listing the penetrance
                             table of the model.


           
b.  Running simpen to produce a model and data:
    
    Command:
    simpen sampleprob.simpen prob

    Input file parameters:
    
    # parameters for model creation
    SEED        1000
    LOCI        2
    FREQ        0.8 0.2
    TABLEVAR    0.10
    USETABLEVAR OFF
    HERIT       0.10
    HERITWEIGHT 5
    MARGVAR     0.00001
    MARGWEIGHT  8
    MODELS      1
    GEN         3000
    POPSIZE     50
    DEMES       5
    MUTATE      0.01  
    CROSS       0.6
    VERBOSE     ON
    SUBMODELS   ON
    UPDATE      100
    SIMTYPE     prob

    # parameters for data simulation
    SIMSETS     1
    AFFECTED    100
    UNAFFECTED  100
    SIMLOCI     20
    ALLELELIMITS 0.05 0.5
    PHENOCOPY   0.0
    GENOERROR   0.0 
    
    Output files:
    prob.out       -->  Lists the penetrance of the models evolved and the 
                        stats.
    prob.1.datasim -->  File used with data simulator.  The '1' refers to
                        this using the first model evolved by simpen.
    prob.1.smod    -->  File used with data simulator listing the penetrance
                        table of the model.
    prob.1.1.txt   -->  Simulated data file.  Lists the affected status
                        followed by the genotypes for every simulated
                        individual.  The '.1.1. extension refers to this
                        being generated by the first model (prob.1.smod)
                        and being the first dataset simulated.  If SIMSET
                        was set to 2 there would be a prob.1.2.txt file also.
    prob.1.loc     -->  Lists the locus information:  allele frequencies and
                        recombination fraction 
  
          
    
c.  Running simpen using the datasim file produced or created by the user.

    Command:
    simpen prob.1.datasim prob
    
    Input file parameters:
    
    # this file is needed for library data simulation config file
    RAND 1000
    MODELFILES 1
    prob.1.smod 1.0
    SIMSETS 1
    GENOTYPEERROR 0
    PHENOCOPY 0
    AFFECTED 100
    UNAFFECTED 100
    SIMLOCI 10
    ALLELELIMITS 0.05 0.5
    DEFAULTALLELE 0.8 0.2
    ALLELEFREQS
    4 0.8 0.2
    9 0.8 0.2
    
    
    Output files:
    prob.1.txt     --> Simulated data file.  Lists the affected status
                       followed by the genotypes for every simulated
                       individual.  The '.1.txt' extension shows this is
                       the first dataset simulated.  If SIMSET was set to 2 
                       there would be a prob.2.txt file also.
    prob.loc       --> Lists the locus information:  allele frequencies and
                       recombination fraction 


3.  Input files

The user must provide a configuration file listing the parameters for use in the 
program.  Parameters should be capitalized as shown below.  Most parameters have
default values as indicated.  All parameters should list the parameter name
followed by one or more spaces and then the parameter values.  Each parameter
should be on its own line.

###########  Simpen file parameters ###########

Parameter     Default    Description
---------     ------    -----------------
SEED          Time       Seed for random generator

LOCI          2          Number of loci in model

FREQ          0.8 0.2    These are the frequencies of the alleles. A 
                         locus only has two alleles in the current 
                         version of simpen.  The user can specify 
                         different allele frequencies for each locus.  
                         If there are more loci than specified 
                         frequencies, all extra loci have the same 
                         frequencies as the first listed.  In the 
                         example all 4 loci would have allele 
                          frequencies set to 0.7 and 0.3. 

TABLEVAR      0.1        Target table variance.

USETABLEVAR   ON         When set to ON, table variance will be used
                         as part of the fitness function.  Use "OFF"
                         to ignore table variance.

HERIT         0.1        Target heritance always used in fitness function
                                    for selection of models.
                                    
HERITWEIGHT   5          Weight given in fitness function to 
                         difference between model's heritability and the
                         target heritability                                    

MARGVAR       0.00001    Marginal variance target always used in fitness
                         function.
                                    
MARGWEIGHT    5          Weight given in fitness function to 
                         difference between model's marginal variance
                         and the target marginal variance.
                                    
PENTARGET     0.0        Marginal penetrance target used in fitness
                         function if this parameter is included
                         in the configuration file.  Otherwise,
                         the marginal penetrance will not be 
                         included in the fitness function.

GEN           5000       Maximum number of generations.  The program
                         will also terminate if a solution has
                         met the fitness criteria.

POPSIZE       100        Number of individuals in each deme

DEMES         10         Number of demes used in genetic algorithm
 
CROSS         0.9        Crossover rate for the genetic algorithm

MUTATE        0.001      Mutation rate per element in the genome.

MODELS        1          Number of models to generate.

SUBMODELS     OFF        When on, the program calculates the
                         statistics for submodels contained within
                         the best model found.  For model ABC,
                         it would calculate the statistics for
                         AB, AC, and BC.

VERBOSE       OFF        Reported output will include additional
                         information such as the marginal penetrances
                         marginal variance, table variance, 
                         heritability, and fitness for each model.
                                    
UPDATE        100        Interval at which the current best model
                         fitness is displayed.

SIMTYPE      prob        Must be either 'pop' or 'prob'.  'pop' 
                         means simulation will be population-
                         based while 'prob' will be probability-
                         based.
                         
SIMSETS         0        When not set or set to zero, no simulated data
                         will be produced.  Only the .out, .datasim, and
                         .loc files will be produced.  If set >= 1, then
                         .txt files are produced with the number of sets
                         specified here for each model.
                         
AFFECTED      200        Number of affected individuals in a probability
                         based dataset.  Ignored for population based.

UNAFFECTED    200        Number of unaffected individuals in a probability
                         based dataset.  Ignored for population based.
                         
SIMLOCI        20        Total number of loci in a simulated probability
                         based dataset.  Only the loci in the simpen model
                         will affect individual disease status.  Not used in 
                         population based simulation.
                         
PHENOCOPY     0.0        Phenocopy rate in the simulated set.  Used in both
                         probability and population based simulation.
                         
GENOERROR     0.0        Genotype error rate in the simulated set.  Used in both
                         probability and population based simulation.
                   
ALLELELIMITS  none       Sets the range of the minor allele frequency in a
                         simulated SNP.  Used in both probability and 
                         population based simulation.  When it is not used
                         all SNPs will have the frequencies specified by the
                         FREQ keyword above.

SIMPOPSIZE     200       Total size of the simulated population in a population
                         based data simulation.  Ignored for probability based.
                         
NUMGENS         10       Number of generations for population based data 
                         simulation to run.  Ignored for probability based.
                         
GENES            5       Number of genes to simulate in population based data
                         data simulation.  Ignored for probability based.
                         
MINSNP           5       Minimum # of SNPs per gene in population based data
                         simulation.  Ignored for probability based.
                         
MAXSNP           5       Maximum # of SNPs per gene in population based data
                         simulation.  Ignored for probability based.
                         
MINRECOMB    0.005       Minimum recombination rate between SNPs in a gene
                         for population based data simulation.  Ignored
                         for probability based.
                         
MAXRECOMB     0.05       Maximum recombination rate between SNPs in a gene
                         for population based data simulation.  Ignored
                         for probability based.



###########  Datasim file parameters ###########

Parameter     Default    Description
---------     -------    -----------------    
RAND             time    Set random seed 

MODELFILES       none   Lists number of model files used.  The path to each
                         model file should follow with one per line.  After
                         the path, the fraction of individuals using this model
                         should be listed.  An example where 2 model files have
                         been used:
                         MODELFILES 2
                         /home/dudek/prob.1.smod 0.7
                         /home/dudek/prob.2.smod 0.3

SIMSETS             1    Number of datasets to create.  Each dataset will be in
                         its own file.
                         
GENOTYPEERROR     0.0    Genotype error rate in simulated data.

PHENOCOPY         0.0    Phenocopy rate in simulated data.

AFFECTED          200    Number of affected individuals in a probability
                         based dataset.
                         
UNAFFECTED        200    Number of affected individuals in a probability
                         based dataset.
                         
SIMLOCI            20    Total number of loci in a simulated probability
                         based dataset.  Only the loci in the simpen model
                         will affect individual disease status.  Not used in 
                         population based simulation.
                         
ALLELELIMITS     none    Sets the range of the minor allele frequency in a
                         simulated SNP.  The frequencies for all loci not 
                         listed in ALLELEFREQS will be set randomly using this
                         range.  Used in both probability and 
                         population based simulation.  When it is not used
                         all SNPs will have the frequencies specified by the
                         DEFAULTALLELE keyword below.                     

DEFAULTALLELE 0.5 0.5    Default allele frequencies.  Used if ALLELELIMITS not
                         set.

ALLELEFREQS      none    Specifies the allele frequencies for individual SNPs
                         in the simulation.  The penetrance model loci should
                         be set here to match the frequencies used in generating
                         the model.

POPSIZE           200    Total size of the simulated population in a population
                         based data simulation.  Ignored for probability based.
                         
NUMGENS            10    Number of generations for population based data 
                         simulation to run.  Ignored for probability based.
                         
GENES               5    Number of genes to simulate in population based data
                         data simulation.  Ignored for probability based.
                         
MINSNP              5    Minimum # of SNPs per gene in population based data
                         simulation.  Ignored for probability based.
                         
MAXSNP              5    Maximum # of SNPs per gene in population based data
                         simulation.  Ignored for probability based.
                         
MINRECOMB       0.005    Minimum recombination rate between SNPs in a gene
                         for population based data simulation.  Ignored
                         for probability based.
                         
MAXRECOMB        0.05    Maximum recombination rate between SNPs in a gene
                         for population based data simulation.  Ignored
                         for probability based.



                  
4. Fitness function 

The genetic algorithm uses a fitness function to evolve a model.  The function
evaluates the fitness of each model in the population.  The fitness is 
dependent on the parameters supplied in the configuration file.  Marginal 
penetrance variance and heritability always affect fitness.  Table variance
and marginal penetrance target only affect when set in configuration file.
Maximum fitness is 1.0.  The genetic algorithm terminates when it finds a model
with fitness = 1.0 or the maximum number of generations is reached.

Marginal Penetrance Variance -- If the marginal penetrance variance is greater
than the MARGVAR parameter in the configuration file, fitness is reduced by the 
weighted difference of the two.

(|Marginal Penetrance Variance - MARGVAR|) * MARGWEIGHT


Heritability -- If the model's heritability is more than .002 different
from the HERIT parameter, then the fitness is reduced by the weighted 
difference of the two when this parameter has been set.

(|Heritability - HERIT|) * HERITWEIGHT


Odds Ratio -- If the model's odds ratio is more than .05 different
from the ODDSRATIO parameter, then the fitness is reduced by the weighted
difference of the two when this parameter has been set.

(|Odds Ratio - ODDSRATIO|) * ODDSWEIGHT


Target marginal penetrance -- If the PENTARGET parameter is set in the 
configuration file and the average marginal penetrance is more than .003
different from the PENTARGET parameter, the fitness is reduced by the
difference of the two.

(|Average Marginal Penetrance - PENTARGET|) 


Table cell variance -- If the USETABLEVAR parameter is set to on and the 
table cell variance is less than TABLEVAR, fitness is reduced by the difference
of the two.

TABLEVAR - table cell variance



5. Genotype data simulation

Simpen can simulate genotype data using probability or population based 
approaches.  In the probability based approach, the user indicates the 
number of affected and unaffected individuals desired.  The program randomly
determines each genotype based on the allele frequencies and compares the 
individual's genotype at the model loci against the penetrance table.  
Affected and unaffected individuals are selected as needed.

In a population based approach an initial population is generated without
regard to affection status.  The program creates a map based on the parameters
specified and the population undergoes a user-specified number of generations.
The program produces each succeeding generation by randomly selecting 
individuals and crossing them.  The individuals undergo crossover and each 
pass a chromosome to the offspring.  Affection status is assigned after
the last generation.



6.  Screen Output

While processing, the program displays progress to standard output.  If VERBOSE
is off, the only number of generations and current fitness of the best model
are shown as below:

 Count        Fitness
---------------------
   [0]      0.0359932
 [100]       0.724547
 [200]       0.831487
 [300]       0.911068
 [400]       0.954538
 [500]              1
 
If VERBOSE is on, standard output shows the marginal penetrance 
variance, table cell variance, heritability, average marginal penetrance, and
fitness as below:

 Count   MarginalVar      TableVar  Heritability   MarginAvg  Odds Ratio   Fitness
---------------------------------------------------------------------------------------------
   [0]      0.011645     0.0880666      0.193145    0.230875     4.92276  0.774428
 [100]     0.0030982     0.0354254      0.125885   0.0999134     4.41225  0.963467
 [200]   0.000438966      0.040634      0.082714    0.102099     5.02483  0.996496
 [300]    9.2785e-05     0.0583515      0.095366    0.101905     4.96148  0.999266
 [400]   9.11248e-05     0.0510216     0.0942105    0.101798      4.9502  0.999279
 [500]   3.54199e-05     0.0717049     0.0966821    0.100679     4.97871  0.999725
 [600]   3.15381e-05     0.0722288     0.0976447     0.10035     5.02085  0.999756
 [700]   2.32348e-05     0.0482382     0.0949537    0.100415     5.03098  0.999822
 [800]   2.13617e-05      0.042171     0.0950449    0.100267     5.04817  0.999837
 [900]   1.85207e-05     0.0403882     0.0941406    0.100156     5.03602   0.99986
[1000]   1.42806e-05     0.0337783     0.0929668    0.100496     4.95299  0.999894
[1100]   1.06053e-05     0.0384778     0.0945335   0.0997367     5.04992  0.999923
[1200]   7.76386e-06     0.0337406     0.0918796   0.0991904     4.99375  0.999946
[1300]   5.95695e-06     0.0343497     0.0926699    0.100123     5.03312   0.99996
[1400]   5.45445e-06      0.056188     0.0934125    0.100013     5.04375  0.999964
[1500]   4.12609e-06     0.0531033     0.0896092   0.0992318     4.95617  0.999975
[1600]   2.07659e-06     0.0554497      0.089766    0.099108      4.9612  0.999991
[1700]   1.80956e-06      0.051187     0.0899207   0.0990485     4.96251  0.999994

 


 
7.  File Output

The program produces multiple output files. The files have the extensions .out, 
.smod, .datasim, .txt and .loc.  The names are dependent on the optional second 
command line argument to the program.  The default names are models.out, models.smod, 
and models.datasim.  The .txt and .loc files are only produced if the program
simulates genotype data.


.out file -- Displays the models produced by the program.  It always displays 
penetrance values for the genotypes of the model as below:

-----------------------------------
Model #1
-----------------------------------
AABBCC 0.084
AABBCc 0.050
AABBcc 0.107
AABbCC 0.055
AABbCc 0.343
AABbcc 0.000
AAbbCC 0.112
AAbbCc 0.078


If VERBOSE is on, it also displays the marginal penetrances and other information
on the model.


Marginals = 0.0977040 0.0998368 0.0992048 0.0974544 0.100419 0.0985392 0.0980720 0.0991648 0.0986928

Model          MarginalVar    TableVar    Heritability  MarginAvg   Odds Ratio  Fitness
---------------------------------------------------------------------------------------
A B C          9.51708e-07    0.0562089   0.0914816     0.0984465   5.03461     1.00000



If SUBMODELS is on, simpen calculates the results for each multilocus submodel.
For example, for ABC it calculates AB, BC, and AC and displays the results.

Model          MarginalVar    TableVar    Heritability  MarginAvg   Odds Ratio  Fitness
---------------------------------------------------------------------------------------
A B C          9.51708e-07    0.0562089   0.0914816     0.0984465   5.03461     1.00000
A B            1.38382e-06    0.0196430   0.0230279     0.0984465   2.44956     0.948988
A C            6.22398e-07    0.0161395   0.0231691     0.0984465   2.42863     0.948573
B C            1.02802e-06    0.00435284  0.0379355     0.0984465   3.71467     0.974293




.smod file -- The smod file is used in the data simulation software.  Every
model simulated produces a smod file.  The number of the model is incorporated
into the file name.


# specify which loci are disease loci
DISEASELOCI 4 9 

# Frequencies used in creating this model: 0.8 0.2
PENTABLE
AABB 0.309
AABb 0.541
AAbb 0.486
AaBB 0.550
AaBb 0.080
Aabb 0.171
aaBB 0.431
aaBb 0.284
aabb 0.517


.datasim file -- Also used by the data simulation software.  The file lists
all the parameters necessary for producing a dataset with the model specified
in the file.  

RAND 1000
MODELFILES 1
prob.1.smod 1.0
SIMSETS 1
GENOTYPEERROR 0
PHENOCOPY 0
AFFECTED 100
UNAFFECTED 100
SIMLOCI 10
ALLELELIMITS 0.05 0.5
DEFAULTALLELE 0.8 0.2
ALLELEFREQS
4 0.8 0.2
9 0.8 0.2


.txt file -- Only produced when data are simulated.  First column is affection
status and each additional column is the genotype of the individual at that
locus.

1 1 1 0 1 1 1 0 0 1 0
1 1 1 0 1 0 0 0 0 0 0
1 0 1 1 1 0 1 1 1 1 1
0 0 1 0 1 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0 0
1 1 1 0 1 0 0 0 1 1 1
0 1 0 1 1 0 0 0 0 0 0
0 0 2 0 1 0 0 1 0 0 1
1 0 2 0 1 0 0 1 0 0 0
1 0 0 0 0 0 0 0 1 0 0


.loc file -- Only produced when data are simulated.  Displays the allele 
frequencies and the recombination fraction for the loci used in the simulation.
The recombination fractions for loci used in a probability based simulation
are zero.

Loc All1 All2 Recomb
  0 0.75 0.25 0.500
  1 0.82 0.18 0.010
  2 0.86 0.14 0.029
  3 0.70 0.30 0.015
  4 0.80 0.20 0.500
  5 0.53 0.47 0.015
  6 0.63 0.37 0.032
  7 0.85 0.15 0.050
  8 0.93 0.07 0.500
  9 0.80 0.20 0.027
                                                                        


