```{r setup, include=FALSE}
knitr::opts_chunk$set(echo = TRUE)
```

# Genetic drift simulation in R

## Step 1 - initial simulation

#### Set up the initial population to begin in Hardy-Weinberg equilibrium.

First, we will need to make a data frame that has the initial parent population of males in it.

```{r make.parent.male.pop}



```

To make the females, just assign the males into an object called females:

```{r make.parent.female.pop}



```

Combine the males and females, and "unlist" them to make it easy to count alleles across both columns.

```{r combine.and.unlist.parents}



```

Now we can count up how many alleles are present in the population:

```{r count.unique.alleles}



```

We can get the heterozygosity for the parent population by comparing the alleles:

```{r heterozygosity}



```

Calculate the allele frequencies for the parent alleles:

```{r allele.frequencies}



```

Set up a blank data frame that we can use to collect the simulation results.

```{r make.results.data.frame}



```

Now put the outputs into a single vector:

```{r bundle.output}



```

Append the row to the results table:

```{r add.initial.to.results}



```

Now we can simulate reproduction. First, randomly select the male and female breeders:

```{r select.breeders}



```

Select the alleles from each parent, and make the male offspring:

```{r make.offspring}



```

Replace the parents with this generation's offspring:

```{r offspring.mature.to.adults}



```


Repeat using a loop:

```{r loop}



```

Make a graph of the allele frequencies:

```{r allele.freq.graph}



```

Make a graph of heterozygosity:

```{r hz.graph}



```

Make a graph of allelic diversity:

```{r allelic.graph}



```

## Second model - repeat the first model 100 times

Make a function out of the previous simulation:

```{r simulation.run.function}



```


Use an lapply() command to run the simulation 100 times:

```{r lapply.repeat.100.times}



```

Use an sapply() command to get the generation an allele went to fixation for each simulation:

```{r sapply.to.get.fixation.gen}



```

Use another sapply() command to get which allele it was that went to fixation:

```{r sapply.to.get.allele.fixed}



```

Combine the gen.fixed and allele.fixed into a data frame:

```{r fixation.data}



```

Get some statistics on the rate of fixation:

```{r fixation.stats}



```

Plot the generation fixed as a boxplot.

```{r gen.fixed.boxplot}



```