---
title: "350.debugging"
author: "William Revelle"
date: "05/04/2024"
output: html_document
---

```{r setup, include=FALSE}
knitr::opts_chunk$set(echo = TRUE)
options(width=100)
```
# An example of function evaluation and testing


The lm function is the workhorse for linear regression. It works on raw data and provides very useful diagnostics.

*lmCor*,  originally named *setCor*, was developed to handle input from correlation matrices. It was eventually elaborated to mimic much of the functionality of *lm*, and to work with raw data as well as correlations, but with various options added.  In the process of doing so, it was necessary to check against *lm* for various options. 

(In 2020, Robert Stupinsky reported an error with zero centering and standarizing. Thus, this check of the revised code was written. It seems to have worked.  See the news file for version 2.0.12.)

We first do a simple multiple regression with three predictors and compare that to the lmCor results. By default, lmCor standardizes and zero centers the data. We turn this off for this first example.

Here we go through a number of tests.  We need to compare the output carefully to detect errors.  

We first use a small data set built into core-R: `attitude`.  We then test it on a larger and more complicated data set: `sat.act`.



# The data

```{r}
library(psych)
library(psychTools)
describe(attitude)
R <- lowerCor(attitude) #find the correlations for further testing
```

# Simple regression

We need to specify some of the options for `lmCor` to make the results match `lm`. In addition, we need to wrap the `lmCor` commands with a print statement to get the same number of decimal places as `lm`

```{r simpleregression}
summary(lm(rating ~ complaints + privileges + learning, data=attitude))
print(lmCor(rating ~ complaints + privileges+ learning, data=attitude,std=FALSE,zero=FALSE),digits=4)
lmCor(rating ~ complaints + privileges+ learning, data=attitude) #default values
lmCor(rating ~ complaints + privileges+ learning, data=R,n.obs=30)
```

# Linear regression with an interaction

To find interactions, we express it with a product.  Unless we zero center, this changes the interpretability of the first order terms.

```{r}
mod1 <- lm(rating ~ complaints + privileges + learning, data=attitude)
mod2 <- lm(rating ~ complaints * privileges + learning, data=attitude)
#now show them
summary(mod1)
summary(mod2)
#zero center
mod3  <- lm(rating ~ complaints * privileges + learning, data=data.frame(scale(attitude,scale=FALSE)))
summary(mod3)
anova(mod2,mod3)
```                        

## Compare to lmCor
```{r}
sc1 <- lmCor(rating ~ complaints + privileges + learning, data=attitude)
sc2 <- lmCor(rating ~ complaints * privileges + learning, data=attitude)
summary(sc1)
summary(sc2)
anova(sc1,sc2)
```