EZM_PNN.PV_Regressions_Rearing.Rmd

---
title: "PNN/PV Regressions Stats/Data + Rearing Interaction"
author: "Alissa Valentine and Kelsea Gildawie"
---
# Loop Notes:
## - run the regressions with every behavioral measure (3) with every PNN/PV measure (4 each for PL and IL)
## - Algorithm: [behavioral measure] ~ [centered PNN/PV measure] + Rearing + Housing + [centered PNN/PV measure]:Rearing (separately for males and females)
## - 48 regressions (3 behavior measure, 4 PNN/PV measures, 2 regions, 2 sexes)

# Loading library, saving dataset, exploring data
```{r}
# Loading library
library(foreign)
library(jtools)
library(interactions)
library(leaps)
library(dplyr)
# Saving SPSS file from my local computer (please update this address depending on where you save this file)
dataset <- read.spss(file = "/Users/alissavalentine/EZM_PNN_Regressions_min2.sav", to.data.frame=TRUE)
# Reformating as a tibble
dataset <- as_tibble(dataset)
# checking out data values
summary(dataset)
```

## Creating Male and Female datasets
### Removing outlier - subject 70.52
```{r}
# Creating data set for males, filtering out outlier in subject 70.52
dataset_males <- dataset %>%
  filter(dataset$'Sex' == "Male", dataset$'Subject' != "70.52")
# Creating data set for females
dataset_females <- dataset %>%
  filter(dataset$'Sex' == "Female")
```

# *Preparing Variables for Running Multiple Regression*

## Creating new variables that hold the names of behavior and PNN/PV variables
```{r}
# Saving column names for the behavior variables into 'behavior' list
behavior <- colnames(dataset_females[5:7])
# Saving column names for the PNN/PV variables into 'PNNPV' list
PNNPV <- colnames(dataset_females[8:15])
```

# *Running the Female Multiple Regression*

```{r}
# Female multiple regression loop using the following equation:
# linear model = *behavior measure*  ~  *PNN/PV measure* + Rearing + Housing + *PNN/PV measure*:Rearing

# This loop iterates through each behaviorial measure, running a mulitple linear regression with behavioral variable as DV, and each PNN/PV variable as the IV being iterated through the regression for each behavior variable (i.e. 3 behaviors x 8 PNN/PV variables = 24 regressions)
for (i in 1:length(behavior)){ # iterate through each behavior
  b <- behavior[i] # saving behavior variable name as 'b' for each iteration
  assign((paste("femalelm", b, sep = "_")), lapply(PNNPV, function(x){
    lm(as.formula(paste(b," ~ ", x, "*Rearing", sep = "")), data = dataset_females)})) # assigning the results from the multiple linear regression (where PNN/PV variables are iterating) to a list named 'femalelm' + name of behavioral measure
}
```

# *Running the Male Multiple Regression*

```{r}
# Male multiple regression loop using the following equation:
# linear model = *behavior measure*  ~  *PNN/PV measure* + Rearing + Housing + *PNN/PV measure*:Rearing

# This loop iterates through each behaviorial measure, running a mulitple linear regression with behavioral variable as DV, and each PNN/PV variable as the IV being iterated through the regression for each behavior variable (i.e. 3 behaviors x 8 PNN/PV variables = 24 regressions)
for (i in 1:length(behavior)){ # iterate through each behavior
  b <- behavior[i] # saving behavior variable name as 'b' for each iteration
  assign((paste("malelm", b, sep = "_")), lapply(PNNPV, function(x){
    lm(as.formula(paste(b," ~ ", x, "*Rearing", sep = "")), data = dataset_males)})) # assigning the results from the multiple linear regression (where PNN/PV variables are iterating) to a list named 'malelm' + name of behavioral measure
}
```

# *Accessing One Linear Model at a Time*

# Key to female/male linear model lists:
## [[1]] - PV_count_PL_c
## [[2]] - PNN_pos_PV_count_PL_c
## [[3]] - PNN_pos_PV_intensity_PL_c
## [[4]] - PV_pos_PNN_intensity_PL_c
## [[5]] - PV_count_IL_c
## [[6]] - PNN_pos_PV_count_IL_c
## [[7]] - PNN_pos_PV_intensity_IL_c
## [[8]] - PV_pos_PNN_intensity_IL_c

# *Statistics for Female Subjects*

# Female - Crossings
```{r}
# Returning summary, summ, and anova stats for every linear model in the femalelm_crossings list
lapply(femalelm_crossings, function(x){
    print(summary(x))
    print(summ(x))
    print(anova(x))
    plot(x)
})
```

# Female - Frequency in the open
```{r}
# Returning summary, summ, and anova stats for every linear model in the femalelm_freq_open list
lapply(femalelm_freq_open, function(x){
    print(summary(x))
    print(summ(x))
    print(anova(x))
    plot(x)
})
```

# Female - Head poke duration
```{r}
# Returning summary, summ, and anova stats for every linear model in the femalelm_head_poke_duration list
lapply(femalelm_head_poke_duration, function(x){
    print(summary(x))
    print(summ(x))
    print(anova(x))
    plot(x)
})
```

# *Statistics for Male Subjects*

# Male - Crossings
```{r}
# Returning summary, summ, and anova stats for every linear model in the malelm_crossings list
lapply(malelm_crossings, function(x){
    print(summary(x))
    print(summ(x))
    print(anova(x))
    plot(x)
})
```

# Male - Frequency in the open
```{r}
# Returning summary, summ, and anova stats for every linear model in the malelm_freq_open list
lapply(malelm_freq_open, function(x){
    print(summary(x))
    print(summ(x))
    print(anova(x))
    plot(x)
})
```

# Male - Head poke duration
```{r}
# Returning summary, summ, and anova stats for every linear model in the malelm_head_poke_duration list
lapply(malelm_head_poke_duration, function(x){
    print(summary(x))
    print(summ(x))
    print(anova(x))
    plot(x)
})
```