This is an [R Markdown](http://rmarkdown.rstudio.com) Notebook. When you execute code within the notebook, the results appear beneath the code. 



```{r}
library(broom)
library(tidyverse)
```


# qlogis is the logit function
# plogis is the logistic function

```{r}
ys <- seq(from = 0, to = 1, by = 0.001)
data <- tibble(ys)
data <- mutate(data, xs = qlogis(ys), zs = plogis(xs) )
tibble(data)
```
```{r}
plot(data$xs, data$ys, pch=19)
```

Let's get the crabs dataset. We will try to predict whether a female horseshoe crab has satellites or not based on its width. Satellites are creepy male crabs that linger nearby.

```{r}
crabs <- read_table('http://users.stat.ufl.edu/~aa/cat/data/Crabs.dat')
#crabs <- mutate(crabs, widthC = width - mean(width), widthZ = scale(width))
tibble(crabs)
```


What does the y ~ width data look like?
```{r}
plot(crabs$width, crabs$y, pch=19)   # (plot(x,y) puts x on x-axis and y on y-axis)
```


```{r}
mdl <- glm(y ~ width, data = crabs, family = 'binomial')
summary(mdl)
```
What do these coefficients mean?
They are in log odds. The width coefficient is positive so the larger the width there is an increase in the log odds of a satellite.
```{r}
intercept <- tidy(mdl)$estimate[1]
slope <- tidy(mdl)$estimate[2]

intercept + slope * 20 # width is 20cm
intercept + slope * 30 # width is 30cm
```

These are the log odds. What does that mean in probabilties?

```{r}
plogis(intercept + slope * 20) # width is 20cm
plogis(intercept + slope * 30) # width is 30cm
```




What do the model predictions look like?
```{r}
widthVals <- seq(from=10, by=0.5, to=40)
mypreds <- tibble(width = widthVals)
mypreds$fit <- predict(mdl, newdata = mypreds)
mypreds
plot(mypreds$width, mypreds$fit, pch=19)
```

So now let's take logistic function to get back a score between 0 and 1.
```{r}
mypreds <- mutate(logistic = plogis(fit), mypreds)
mypreds
plot(mypreds$width, mypreds$logistic, pch=19)
```

Let's plot the scatterplot and model together using ggplot2.

```{r}
ggplot(crabs, aes(x = width, y = y)) +
geom_point(size = 2, alpha = 1) +
geom_line(data = mypreds, aes(x = widthVals, y = logistic)) +
theme_minimal()
```




```{r}
ggplot(crabs, aes(x=width, y=y)) +
  geom_point() +
  geom_smooth(data = crabs, aes(x = width, y = y),
    method = "glm", method.args = list(family = "binomial"), se=FALSE)
```

