```{r setup, include=FALSE}
knitr::opts_chunk$set(echo = TRUE)
download.file('https://wkristan.github.io/template.docx', 'template.docx', mode = 'wb')
```

## Screening data for use in multivariate analysis

First, import the shells data.

```{r import.shells}



```

Check if the file has complete data (no missing cells):

```{r complete.cases.shells}



```

## Initial graphical screening

We will need to make boxplots, normal probability plots, and histograms for each variable to check for outliers and distributional issues.

Since we are screening multiple variables we will lay out the plots in a single multiple-panel graph so that we can more easily look at all of the variables at once. Set up the graphical parameters (par()) for 2 rows of graphs with 3 columns, and then plot the boxplots (don't worry about the results='hide' in the curly braces, this just suppresses text output from the graph functions you use, but the graphs are still displayed).

```{r shell.boxplots, results='hide'}



```

#### Question: how can you detect a skewed distribution in a boxplot?

> 

For the normal probability plots and histograms we will need to calculate residuals, which will be easier to do if we make a vector that has the names of the numeric variables in it:

```{r make.shell.var.list}



```

We will want to use residuals for normality plots to get rid of differences between group means. We will get those by fitting a MANOVA first, and then extracting residuals from the fitted model.

```{r residuals.from.MANOVA}



```

Now make normal probability plots with the residuals:

```{r shell.qqplot, results='hide'}



```

#### Question: how do you know from a normal probability plot if the variable is normally distributed?

> 

#### Question: which variable is the most right-skewed, according to your boxplots and normal probability plots?

> 

We can check the distribution of the residuals using histograms as well:

```{r shells.histo, results='hide'}



```

#### Question: is the same variable you thought was right-skewed based on boxplot and normal probability plot also right skewed as a histogram? How can you tell skew from a histogram?

> 

Next we can make a scatterplot matrix that puts scatterplots in the lower triangle, and correlations in the upper triangle. Run this code chunk to add the panel.cor() function into your environment:

```{r panel.cor.function}

panel.cor <- function(x, y, digits = 2, prefix = "", cex.cor, ...)
{
    par(usr = c(0, 1, 0, 1))
    r <- abs(cor(x, y))
    txt <- format(c(r, 0.123456789), digits = digits)[1]
    txt <- paste0(prefix, txt)
    if(missing(cex.cor)) cex.cor <- 0.8/strwidth(txt)
    text(0.5, 0.5, txt, cex = cex.cor * r)
}

```

We can make the scatterplot matrix using the pairs() function:

```{r pairs.plot}



```

#### Question: which variable is involved in the plots that show non-linear relationships?

> 

Weight seems to be right-skewed and non-linearly related to other variables, so let's try log-transforming it:

```{r log.transform.weight}



```

Add the residuals for ln.weight to shells.resid:

```{r log.transformed.residuals}



```


Now make a boxplot and qqPlot for ln.weight.

```{r replot.ln.weight, results='hide'}



```

#### Question: did log-transforming weight help with the distribution?

> 

Make a list of the variables you will use in the pairs() plot, with ln.weight taking the place of weight:

```{r shell.variables}



```

Plot the scatterplot matrix using log.weight instead of weight.

```{r pairs.plot.log.weight}



```

#### Question: did the transformation straighten out the curved relationships with weight? Did it introduce any new issues that weren't there before?

> 

## Multivariate distribution check

Now we will get the residuals from a MANOVA of the shells data, and check if their Mahalanobis distances follow the Chi-square distribution, like we expect them to.

Calculate the means of the residuals (these should all be 0 or very close to it):

```{r means.of.residuals}



```

The covariance matrix of the residuals next:

```{r covariance.matrix}



```

Mahalanobis distances are calculated from the residuals, the vector of means, and the covariance matrix:

```{r mahalanobis.distances}



```

Compare the Mahalanobis distances to the Chi-square distribution:

```{r mahalanobis.vs.chisquare, results='hide'}



```

#### Question: based on the graph, are there multivariate outliers?

> 

Make a parallel coordinates plot to help identify which variables are responsible for those big Mahalanobis distance values:

```{r ggparcoord}



```

#### Question: for the observations that are out of the ordinary, which variables are responsible? Look for lines that are in the middle of the distribution for most of the variables, but have one or more measurements that are really high or low.

>

#### Question: we still have the shells, what could we do to fix the problem? If we didn't have the shells anymore what would our options be?

>

## Teeth

Repeat the steps on the buffalo teeth.

Import the data:

```{r import.teeth}



```

Plot boxplots and normal probability plots - note that there is no grouping variable for the teeth, so you don't need to split the variables by groups for the boxplots this time. Put the variable names into an object called teeth.names, and then use lapply() to make the plots:

```{r boxplots.teeth, results='hide'}



```

Now use lapply() to get normal probability plots:

```{r qqPlots.teeth, results='hide'}



```

Now use an lapply() to get histograms:

```{r teeth.histograms, results='hide'}



```

Scatterplot matrix (no need to color by groups this time):

```{r pairs.plot.teeth}



```

### Log transformation - only if needed, skip if not needed

If there are any right-skewed variables log-transform them:

```{r log.transform.teeth}



```

Re-do the boxplots and normal probability plots using the variables you log-transformed in place of their un-transformed versions:

```{r boxplots.normal.probability.plots.log.teeth}


```

Re-do the scatterplot matrix using the log-transformed variables:

```{r scatterplot.matrix.log.teeth}


```

### Mahalanobis distances for teeth data

Calculate Mahalanobis distances - since there are no groups in the teeth data you just need to use colMeans on the variables, and var() to get the covariance matrix for the variables (no need to do a MANOVA and extract residuals):

```{r colMeans.and.var.teeth}



```

Calculate the Mahalanobis distances and add them as a variable in your teeth data:

```{r mahalanobis.distances.to.teeth.dataset}



```

Do a qqPlot on the Mahalanobis distances (remember, df is equal to the number of variables):

```{r mahalanobis.quantile.plot.teeth}



```

#### Question: are there multivariate outliers? Could you see them in the other plots, or was the Mahalanobis distance needed to identify them?

> 

Make a parallel coordinates plot for the teeth data:

```{r teeth.parallel}



```

#### Question: which variables are producing unusual measurements in the teeth data?

>
