# Part 1: Matrix and vector operations. 

1. Solve the following system: 

$$
\left[\begin{array}{ccccc}
a_1 & b_1 &  & & 0 \\
c_1 & a_2 & b_2 & & \\
& \ddots & \ddots & \ddots & \\
& &  & a_{99} & b_{99} \\
0 & & & c_{99} & a_{100} 
\end{array}\right] 
\left[\begin{array}{c}
x_{1} \\
x_{2} \\
\vdots \\
x_{100}
\end{array}\right]
=
\left[\begin{array}{c}
d_{1} \\
d_{2} \\
\vdots \\
d_{100}
\end{array}\right]
$$
where 
$$
a_j = j, \quad b_j = 1/j, \quad c_j = 1, \quad d_j = \sin(j \pi / 200)
$$
and print $x_1, x_2, \dots, x_5$. 

**Solution**:
```{r}
# Define A.
A <- matrix(rep(0, 100*100), nrow = 100, ncol = 100, byrow = TRUE)

for (i in c(1:100)) {
    A[i, i] <- i

    if (i + 1 < 101)  A[i, i + 1] <- 1/i
    if (i - 1 > 0)    A[i, i - 1] <- 1
}

# Define D.
d <- c(1:100)
d <- sin(d*pi/200)

# Solve Ax = d.
x <- solve(A, d)
x[1:5]
```

# Part 2: For loops. 

1. Write a function that uses a `for` loop to calculate the following with a sequence of $m$, and generate a plot for $m$ verses $E_m$. Avoid using a `for` loop, can you complete the same task?  

$$
E_m = 1 + \frac{1}{2} + \dots + \frac{1}{2^m} - \log(2^m)
$$

**Solution**:
```{r}
# Using for loop.
E_m <- function(m) {
  res <- 1
  for (i in c(1:m)) {
    res <- res + 1/(2^i)
  }
  res <- res - log(2^m)
  
  return(res)
}

# Avoid using for loop.
E_m2 <- function(m) {
  res <- 1
  
  index <- c(1:m)
  denom <- 2^index
  res <- res + sum(1/denom)

  res <- res - log(2^m)
  
  return(res)
}

E_m(100)
E_m2(100)
```



# Part 3: Analyze NYC flight delays.   

Install the "nycflights13" package. The data comes from the US Bureau of Transportation Statistics. Using the data, complete the following tasks: 

1. Find all flights that had an arrival delay of >4 hours, return the first 5 row. (Note: `arr_delay` is in mins)
2. Find all flight names that flew from JFK to IAH, i.e. return only unique values of "flight" variable after filtering. Hint: `unique()` would help. 
3. Find how many flights were operated by UA. 
4. Find how many unique flights were operated by UA. 
5. Sort flights that have the most delayed flights. Show the first 5 row. 
6. Generate a scatter plot with x-axis `dist` and y-axis `delay`, where each dot is a unique flights and destination, `dist` is the average distance of each destination `dest`, and `delay` is the average delay time `arr_delay`, with the size of `dot` equals to the count of delay records. 

```{r, warning=FALSE, message=FALSE}
library(tidyverse)
library(nycflights13)
head(flights)
```

**Solution**

1. Find all flights that had an arrival delay of >4 hours, i.e. return the first 5 row. (Note: `arr_delay` is in mins)

```{r}
flights %>% filter(arr_delay > 240) %>% head(5)
```

2. Find all flight names that flew from JFK to IAH, i.e. return only unique values of "flight" variable after filtering. Hint: `unique()` would help. 

```{r}
df <- flights %>% filter(origin == "JFK" & dest == "IAH")
unique(df$flight)
```

3. Find how many flights were operated by UA. 

```{r}
nrow(filter(flights, carrier %in% c("UA")))
```

4. Find how many unique flights were operated by UA. 

```{r}
df <- filter(flights, carrier %in% c("UA")) 
length(unique(df$flight))
```

5. Sort flights that have the most delayed flights. Show the first 5 row. 

```{r}
flights %>% arrange(desc(dep_delay)) %>% head(5)
```

6. Generate a scatter plot with x-axis `dist` and y-axis `delay`, where each dot is a unique flights and destination, `dist` is the average distance of each destination `dest`, and `delay` is the average delay time `arr_delay`, with the size of `dot` equals to the count of delay records. 

```{r}
flights %>%
  group_by(flight, dest) %>%
  summarise(delay = mean(arr_delay), dist = mean(distance), n = n()) %>%
  ggplot() +
  geom_point(aes(x = dist, y = delay, size = n)) 
```


# 5. NYC bus delays

The data is from the kaggle dataset "Bus Breakdown and Delays NYC - When and why a bus was delayed? Bus delays 2015 to 2017". (https://www.kaggle.com/anthobau/busbreakdownanddelays).

The Bus Breakdown and Delay system collects information from school bus vendors operating out in the field in real time. Bus staff that encounter delays during the route are instructed to radio the dispatcher at the bus vendor’s central office. The bus vendor staff are then instructed to log into the Bus Breakdown and Delay system to record the event and notify OPT. OPT customer service agents use this system to inform parents who call with questions regarding bus service. The Bus Breakdown and Delay system is publicly accessible and contains real time updates. All information in the system is entered by school bus vendor staff.

You can find data for years 2015 to 2017.

```{r}
bus.delay <- read_csv("Bus_Breakdown_and_Delays.csv")
```

Alternatively, you can install the `RKaggle` package to load the dataset directly into R.

```{r message=FALSE, warning=FALSE}
# Run this line to install RKaggle
# install.packages("RKaggle")
library(RKaggle)
bus.delay<-get_dataset("anthobau/busbreakdownanddelays")
```

```{r}
head(bus.delay)
```

The following code creates new features 

```{r}
library(lubridate)

bus.delay$occur_date <- as.POSIXct(bus.delay$Occurred_On, format= "%m/%d/%Y %H:%M:%S %p")

# Day of the week
bus.delay$occur_weekday <- wday(bus.delay$occur_date, label = T)
bus.delay$occur_month <- month(bus.delay$occur_date, label = T)
bus.delay$occur_year <- year(bus.delay$occur_date)
#head(bus.delay)
```

1. Group by `occur_month`, and generate a bar plot that showing the number of delays by month. 

```{r message=FALSE, warning=FALSE}
library(tidyverse)
```

```{r}
#library(ggplot2)
# tidyverse
bus.delay %>%
  group_by(occur_month) %>%
  summarise(count = n()) %>%
  ggplot(aes(x = occur_month, y = count)) + #axis
  geom_bar(stat = 'identity') # stacked 
```

The following code generates a new variable called "duration" in mins. 

```{r}
bus.delay$duration <- as.numeric(gsub("([0-9]{1,2}).*$", "\\1", bus.delay$How_Long_Delayed))
```


2. Plot the boxplot of duration. 

```{r}
bus.delay %>%
  ggplot(aes(y = duration)) + 
  geom_boxplot()
```

3. Plot a histogram of durations and group by "Boro".

```{r}
bus.delay %>%
  ggplot() + 
  geom_histogram(aes(x = duration, 
                     fill = Boro), 
                 position = 'dodge', 
                 bins = 3)  
```

4. Plot the density of duration, group by "occur_weekday" and use facet to stratify on" "Boro".

```{r}
bus.delay %>%
  ggplot() +
  geom_density(aes(x = duration,
                   color = occur_weekday,
                   fill = occur_weekday),
               bw = 0.08) +
  facet_wrap(~Boro, ncol = 4)
```

5. Print top five Boro by mean delay.

```{r}
bus.delay %>% 
  group_by(Boro) %>% 
  summarise(mean_delay = mean(duration, 
                              na.rm = TRUE)) %>%
  arrange(desc(mean_delay)) %>%
  head(5)
```