# Purpose

In today's lab, we will discuss how to use R to visualize correlations, calculate correlation coefficients, and assess statistical significance of correlations. You can download the Rmd file [here](lab-1.Rmd) to follow along.

For today's lab, you will need to load the following libraries. 

If you are missing any of the following packages, you can use the following code to install it.

```{r, eval = FALSE}
install.packages("pwr")
```

```{r, message=FALSE, warning = FALSE}
library(tidyverse) # for plotting and data wrangling
library(rio) # for importing data
library(psych) # for covariance and correlation functions
library(apaTables) # for correlation tables
library(pwr) # for power calculations
library(here) # for file paths
library(corrplot) # for heat maps
```


*** 

# Visualizing correlations {#plot}

It is very important to *always* visualize your data. There might be patterns in your data that are not apparent just from looking at a correlation between two variables. We will go over just a few examples here, but there are many different options for visualizing correlations. See [here](http://www.sthda.com/english/wiki/correlation-matrix-a-quick-start-guide-to-analyze-format-and-visualize-a-correlation-matrix-using-r-software#use-corrplot-function-draw-a-correlogram){target="_blank"} if you want to explore more options. 

## Scatter plots

Scatter plots allow us to visualize the relationship between two variables. Here is the code for a scatter plot using `ggplot()`. You can add a line of best fit by adding a `geom_smooth()` layer to the plot. 

```{r, message = FALSE}
ggplot(data = mtcars, aes(x = mpg, y = hp)) +
  geom_point() +
  geom_smooth(method = "lm")

```


## SPLOM plots

"SPLOM" stands for scatter plot matrix. The `pairs.panel()` function from the `{psych}` package allows a quick way to visualize relationships among all the continuous variables in your data frame. The lower diagonal contains scatter plots showing bivariate relationships between pairs of variables, and the upper diagonal contains the corresponding correlation coefficients. Histograms for each variable are shown along the diagonal. 

```{r}
pairs.panels(mtcars, lm = TRUE)
# stars indicates significance
# cex change the size of fonts for the correlation coefficients.
```


## Heat maps

Heat maps are a great way to get a high-level visualization of a correlation matrix. They are particularly useful for visualizing the number of "clusters" in your data if that's something you're looking for. `Thurstone`, built into the `{psych}` package, is a correlation matrix of variables assessing different aspects of cognitive ability. We can plot a heatmap of this correlation matrix using the `corrplot()` function from the `{corrplot}` package. Note: make sure that you are feeding the function a correlation matrix (not the data set).

```{r}
corrplot(corr = Thurstone, method = "square")

```

>**Question:** What do you notice about the structure of this data?

## APA Tables

* The package `{apaTables}` has a very useful function `apa.cor.table()` that creates nicely formatted tables of correlation matrices in APA format. This code prints the table to a word document called "mtcars.doc" that I am storing in the `lab-1` folder of my working directory. `table.number = 1` indicates "Table 1."


```{r, eval = FALSE}
apa.cor.table(mtcars, 
              filename = here("labs", "lab-1", "mtcars.doc"), 
              table.number = 1)
```

<center>
![](images/cars_table.png)
</center>


***

# Covariances and correlations{#corcoeff}

## Covariance

Covariance captures how the variances of two variables are related, i.e., how they *co*-vary. If higher values of one variable tend to correspond with higher values of the other variable, and lower values of one variable tend to correspond with lower values of the other variable, then the covariance would be positive. However, if the two variables are inversely related (i.e., higher values on one variable correspond with lower values on the other variable), then the covariance would be negative.  

$$\large cov_{xy} = {\frac{\sum{(x-\bar{x})(y-\bar{y})}}{N-1}}$$

To calculate covariance, use the function `cov()` from the `{stats}` package. The `cov()` function takes two arguments: the first variable "x" and the second variable "y". 

* Calculate the covariance between `mtcars$mpg` and `mtcars$hp`:

```{r}
cov(x = mtcars$mpg, y = mtcars$hp)
```

* Feeding `cov()` a data frame (of numeric variables) will generate a covariance matrix. Calculate a covariance matrix for the `mtcars` data frame. Round to two decimal places. 


```{r}
cov(mtcars)
```


>**Question:** In the above output, what do the numbers along the diagonal represent? 


## Correlation      

* Correlations are *standardized* covariances. Because correlations are in standardized units, we can compare them across scales of measurements and across studies. Recall that mathematically, a correlation is the covariance divided by the product of the standard deviations of each variable. 

$$\large r_{xy} = {\frac{cov(X,Y)}{\hat\sigma_{x}\hat\sigma_{y}}}$$

* Calculate the correlation coefficient for `mtcars$mpg` and `mtcars$hp` using the `cor()` function (again, from the `{stats}` package).

```{r}
cor(mtcars$mpg, mtcars$hp)
```


* As with covariances, we can generate a matrix of correlations by feeding a data frame to `cor()`. Calculate a correlation matrix for the variables in the `mtcars` data set. Round to two decimal places. 

```{r}
cor(mtcars)
```


* If you are given a covariance matrix, you can convert it to a correlation matrix using `cov2cor()`.

```{r}
# covariance matrix
cov_mat <- cov(mtcars)

# convert to correlation matrix
cov2cor(cov_mat) %>% 
  round(2)
```


***

# Hypothesis testing with correlations{#stats}

Once you have a correlation coefficient, you may want to assess whether the correlation is statistically meaningful. 

* For this example, we are going to work with a data set (n = 60; 30 men and 30 women) that contains variables measuring self-reported conscientiousness (from the BFI) and self-reported physical health. We want to know if there is a significant correlation between conscientiousness and physical health, and then whether or not that correlation differs between men and women.

* Import the data using the following code and check the structure of the data. 

```{r lab-1-25, results='hide'}
health <- import("https://raw.githubusercontent.com/uopsych/psy612/master/labs/lab-1/data/consc_health.csv")

str(health)
head(health)
```

## Visualize the data

```{r, message=FALSE}
ggplot(data = health, aes(x = consc, y = sr_health)) +
  geom_point() + 
  geom_smooth(method = "lm") +
  labs(x = "Conscientiousness", y = "Self-reported health", title = "Conscientiousness and self-reported health")
```

## Statistical hypotheses

$$\large H_{0}: \rho_{xy} = 0$$

$$\large H_{A}: \rho_{xy} \neq 0$$

## Correlation test

You can run a correlation test using the `corr.test()` function from the `{psych}` package.

```{r}
(r_consc_health <- health %>% 
  select(consc, sr_health) %>% 
  corr.test())

```

>**Question:** What can we conclude about the relationship between conscientiousness and health from this data? 

* It is often very useful to save the output of a statistical test to an object that you can then pull useful information out of it. 

```{r}
r_consc_health$ci #confidence interval
r_consc_health$p #p-value matrix
```

## Comparing correlations

We saw above that there is a significant positive relationship between conscientiousness and self-reported health. However, now we want to know whether or not this correlation is significantly different for men and women. 

### Visualizing by group

* First let's plot our data again, but this time split the data by gender. 

```{r}
health %>% 
  ggplot(aes(consc, sr_health, color = gender)) + 
  geom_point(size = 2) +
  geom_smooth(method = "lm") + 
  labs(x = "Conscientiousness", y = "Self-reported Health")
```

>**Question:** What can we conclude from this graph?

### `psych::r.test()`

* To statistically compare the correlations between conscientiousness and self-reported health for men and women, we can use the `r.test()` function. This particular function requires the sample size of the first group and two correlations we're testing against each other. So, we'll need to first run the correlation separately for men and women.

```{r}
health_women <- health %>% 
  filter(gender == "female") 

health_men <- health %>% 
  filter(gender == "male")

r_women <- cor(health_women$consc, health_women$sr_health)

r_men <- cor(health_men$consc, health_men$sr_health)
```

`r.test()` has three arguments:

1. `n` = sample size for group 1 (the female group)
2. `r12` = r for group 1 (the female group) 
3.`r34` = r for group 2 (the male group)


```{r}
r.test(n = nrow(health_women), r12 = r_women, r34 = r_men)
```

>**Question:** What does this test suggest about the correlation between conscientiousness and health across genders?



***

# Minihacks

For these minihacks, you should practice using inline code to write up your answers. Inline code allows you to report a variable outside of a code chunk, e.g., `r round(r_men,2)`.

## Minihack 1

Use the `pwr.r.test` function to answer the following questions: 

a. How many participants should we collect if we expect a small correlation (r = .2) and want to achieve power = .8? Assume alpha = .05.

```{r}
#pwr.r.test(n = , sig.level = , power = , r = )
```

b. How does the number of participants needed to obtain .8 power change when you decrease the significance level to .01? 
```{r}
#pwr.r.test(n = , sig.level = , power = , r = )
```

c. How does the number of participants needed to obtain .8 power change when you increase the effect size (r) to .5? 

```{r}
#pwr.r.test(n = , sig.level = , power = , r = )
```

## Minihack 2

For this minihack, calculate a 95% CI "by hand" (i.e., you can use R but not the `cor.test()` function) for the correlation between conscientiousness and self-reported health from the `health()` dataset. Make sure you get the same answer that was given by `cor.test()`. Hint: when calculating CI's, use the functions `fisherz()` and `fisherz2r()` from `{psych}`.

To review the steps of calculating confidence intervals using Fisher's transformation, see [here](https://uopsych.github.io/psy612/lectures/01-correlation.html#/35){target="_blank"}.

```{r lab-1-45}
# Your code here
```


## Minihack 3

The following table is from [this paper](https://journals.sagepub.com/doi/pdf/10.1177/0956797615586560?casa_token=UgWG-wLWieQAAAAA%3AfOBMH6ehSkk6qp7WtJrABcxdeefuzRJzHKOl6revQ7-d8cT7IfE1ejWPRbUQWmivgaWjp7OyR3LB&){target="_blank"} from Dawtry et al. (2015). 

<center>
![](images/dawtry_table1.png)
</center>

You can import the data using the following code:

```{r lab-1-47}
dawtry_clean <- import("https://raw.githubusercontent.com/uopsych/psy612/master/labs/lab-1/data/dawtry_2015_clean.csv")
```


a. Replicate the correlation matrix from the table using the `apa.cor.table()` function and open it in Microsoft Word. 

```{r}
#Your text here
```

b. Visualize the correlation matrix using a SPLOM plot and/or a heat map.

```{r}
#Your text here
```

