data-analysis.Rmd

---
title: Traffic Data Analysis
output: html_document
---

## Attaching the required libraries
```{r message=FALSE}
library(dplyr)
#options(dplyr.width = Inf)
library(ggplot2)
library(corrplot)
library(knitr)
library(tidyr)
```

## First acquaintance with the data
Loading the data:
```{r cache=TRUE}
df <- read.csv('all-semi-unique.csv')
```
Getting a glimpse of the data to know its dimensions, the data type for each column and some of its assigned values:
```{r}
glimpse(df)
```
Getting the number of NAs for each column:
```{r}
sapply(df, function(x) sum(is.na(x)))
```

# Cleaning the data - Part I (Blind Cleaning)
Removing duplicate rows as well as columns that have a constant value in all observations:
```{r cache=TRUE}
df <- unique(df)
df <- df[sapply(df, function(x) length(unique(x))) > 1]
glimpse(df)
```
There were no duplicate rows but we have successfully reduced the number of columns from 34 to 19. Before being able to clean the data further, we need to make sure that we properly understand it.

# Cleaning the data - Part II (Educated Cleaning)

## `rd.rp.cmid` uniqueness
An initial guess for the meaning of this column is the report's comment ID. We take a look into the uniqueness of this ID.
```{r}
df[df$rd.rp.cmid=="9424155",] %>% nrow
```

We discover that a comment id can be repeated more than once. The reason for this is attributed to different crawling times. It is possible that some data is duplicated if it still appears on bey2ollak.com feed in consecutive crawls. We note that the values for `rd.rp.hr` and `rd.rp.mn` change with respect to the crawl time. This change indicates that these two columns reperesent the elapsed time for this report (with respect to the crawl time of course).  
Such interpretation agrees with bey2ollak.com website since only elapsed times for roads and reports are provided, not the actual times.
This understanding provides two very useful applications:

### Filter rows based on `rd.rp.cmid`
Now, the number of rows can be significantly reduced to represent each report only once   as follows:
```{r}
df <- df[!duplicated(df[,c("rd.rp.cmid")]),]
dim(df)
```
The number of rows is now significantly reduced to 148367.


### Extract report and road absolute timings
Based on our previous understanding of `rd.rp.hr` and `rd.rp.min`, and by observing the website's response data, `rd.hr` and `rd.mn` indicate the elapsed time of the latest report submitted for this report.
So, now we are going to properly represent the time values as time objects indicating the proper timestamp not relative times.
```{r}
df <- df %>% mutate(crawl_date=as.POSIXct(crawl_date, format = "%a %b %e %H:%M:%S UTC %Y", tz = "UTC"))
df <- df %>% mutate(rd.time = crawl_date - rd.mn * 60 - rd.hr*60*60, rd.rp.time = crawl_date - rd.rp.mn * 60 - rd.rp.hr*60*60)
```
Then, we remove the original elapsed time columns:
```{r}
df$crawl_date <- NULL
df$rd.hr <- NULL
df$rd.mn <- NULL
df$rd.rp.hr <- NULL
df$rd.rp.mn <- NULL
```
Finally, we divide into weekday and hour of the day; these two values are the only ones required in our subsequent studies.
```{r}
df <- df %>% mutate(rd.rp.hour=format(as.POSIXct(rd.rp.time, format="%a %b %e %H:%M:%S UTC %Y"), format="%H")) %>% mutate(rd.rp.hour = ((type.convert(rd.rp.hour) + 2) %% 24))
df <- df %>% mutate(rd.rp.wday=format(as.POSIXct(rd.rp.time, format="%a %b %e %H:%M:%S UTC %Y"), format="%a"))
df$rd.rp.time <- NULL
```

## Report status `rd.rp.stid`
It appears from bey2ollak.com response data that the `rd.rp.stid` value indicates the status id of the report being submitted. Values ranging from 1 to 5 correspond to reports on the road congestion, with 1 being the best value (its default `rd.rp.cm` is "7alawa") and 5 the worst value (its default `rd.rp.cm` is "mafeesh amal"). Values from 6 to 10 represent inquiries or special reports like accidents.  
Also, a number of NA values were found in the `rd.rp.stid` column. While a recognizable amount of the comments corresponding to them are either enquiring or acknowledging the existence of radars, other random comments were also found. By examining the ratio of NA values:
```{r}
sum(is.na(df$rd.rp.stid)) / nrow(df)
```
Since they represent only 3% of our data, and their comments are a mix of the other categories, a decision was made to neglect their rows.
```{r}
df <- df[!is.na(df$rd.rp.stid),]
```

## Image columns
By testing bey2ollak.com response data, it was found that values in `rd.rp.img` represent the id of the reporter's profile photo.  
On the other hand, the column `rd.rp.rpImg` represents the id of the attached road image in the report (if any).  
`rd.img` is a boolean value that is absent from most roads in bey2ollak.com response data, which is reflected in its high NA ratio in our dataset:
```{r}
sum(is.na(df$rd.img)) / nrow(df)
```
We choose to remove these column since they will not be useful in our upcoming work.
```{r}
df$rd.rp.img <- NULL
df$rd.rp.rpImg <- NULL
df$rd.img <- NULL
```

## Boolean columns
An interesting observation for the `rd.strq` value is that it is unset only for long travel roads, namely between Cairo & Alex, Cairo & Hurghada, Cairo & Sharm. On the other hand, the majority of the group of roads having this value set are roads inside a city, while a small proportion is shorter travel roads like Zera3y road between Cairo & Banha, Cairo & Tanta, Alex & Tanta. So, there are good reasons to believe that this value is thresholded by the distance. This observation was detected using this snippet:
```{r eval=FALSE}
strqUnset <- df[df$rd.strq==0,c("rd.strq","rd.nm")]
strqUnset <- strqUnset[!duplicated(strqUnset[,c("rd.nm")]),] %>% as.data.frame
strqUnset %>% head
strqSet <- df[df$rd.strq==1,c("rd.strq","rd.nm")]
strqSet <- strqSet[!duplicated(strqSet[,c("rd.nm")]),] %>% as.data.frame
strqSet %>% head
```
No meaningful observations were found for the other two boolean values `rd.new` and `rd.cmrq`, except that most of their values are unset:
```{r}
sum(df$rd.new==0) / nrow(df)
sum(df$rd.cmrq==0) / nrow(df)
```
Since they do not have a meaningful use within our work, and are seemingly unrelated to our metrics of interest, we choose to exclude them:
```{r}
df$rd.strq <- NULL
df$rd.new <- NULL
df$rd.cmrq <- NULL
```

## Final cleaning
Lastly, the values in the `rd.stid` and `rd.time` columns will be unused in our following analyses so we opt to remove them.
```{r}
df$rd.time <- NULL
df$rd.stid <- NULL
glimpse(df)
```

# Interpretation of different column meanings

## Meaning of `rd.ri`
By looking at the html source of bey2ollak.com we can confirm that the each value in `rd.ri` represents an id of the corresponding road. It is interesting to note that some anamolies were detected when running the following snippet.
```{r}
tmpDF <- df %>% group_by(rd.nm) %>% summarize(numRdId=length(unique(rd.ri)))
tmpDF[tmpDF$numRdId > 1,] %>% kable
```

It was found that the two roads displayed in the above results have more than one id, which is weird. We found the reason to be that since bey2ollak.com supports roads in both Cairo and Alex, there are two roads named "Other Roads" in both cities. Also, there are two entries among the two cities for the Sa7rawy road connecting them.

## Reporter name
The two columns `rd.rp.nm` and `rd.rp.fullnm` represent the name and full name of the user who submitted the report, respectively. Two specail names in the `rd.rp.nm` column stand out in term of repititions; **bey2ollakgps** and **fa3el kheir**. The first is some type of an automated agent which gives periodic speed updates, while **fa3el kheir** is the default name for any reporter not willing to supply a specific name.

# Descriptive statistics

## General Summary
Below is a general summary of the data. Here we temporarily convert the road id, report status id, and weekday columns into *factors*, since their initial representations (*int* or *character*) did not allow for any meaningful statistics.
```{r}
tmpDf <- df
tmpDf$rd.rp.cmid <- NULL
tmpDf$rd.ri <- as.factor(df$rd.ri)
tmpDf$rd.rp.stid <- as.factor(df$rd.rp.stid)
tmpDf$rd.rp.wday <- as.factor(df$rd.rp.wday)
tmpDf %>% summary %>% kable
rm(tmpDf)
```

## Type of reports
As previously explained, reports submitted by users can have one of 10 types, as categorized by the `rd.rp.stid` column. While IDs from 1 to 5 describe general info about the current status of the road in terms of congestion, IDs from 6 to 10 are dedicated for special reports. We will now present a general intuition about the type of all the reports submitted by the users, without any distribution across time or space.
```{r}
ggplot(df, aes(factor(rd.rp.stid))) + geom_bar()
```

The following observations are spotted in the bar chart: 

1. The prevailing type is reports describing the 2nd best traffic flow status (whole default value is **lazeez**).
2. A close runner-up is reports of type 10. By running the following investigation on those reports:
```{r}
c1 <- df[df$rd.rp.stid==10 & df$rd.rp.nm=="bey2ollakgps",c("rd.ri","rd.rp.cm")] %>% nrow()
c2 <- df[df$rd.rp.stid==10,c("rd.ri","rd.rp.cm")] %>% nrow()
c1/c2
```
We find that this high number is because the majority of this kind of reports comes from the automated speed update agent.

3. The only other reports having a significant percentage in the dataset are the remaining traffic flow status reports (with *IDs* 1, 3, 4 & 5) and the questions (*ID 10*)

##Report type across time of day
Below is a distribution of each report type across different hours of the day. Each graph represents a report type, ranging from IDs 1 to 10. The x-axis represents the day hour ordered from midnight till 23.00 hr. The y-axis represents the report count.
```{r}
ggplot(df, aes(rd.rp.hour)) + geom_bar() + facet_wrap(~ rd.rp.stid)
```

## Congestion hours
A more useful distribution across day hours should focus only on the reports concerning the traffic flow ( *`rd.rp.stid`=1~5*), neglecting the special reports and questions. Now, the congestion can be represented by a simple mean value across the IDs. The higher this value, the more congested the road is.
```{r}
congestionDf <- df[df$rd.rp.stid < 6,] %>% group_by(rd.rp.hour) %>% summarize(avRep=mean(rd.rp.stid), numRep=n())
ggplot(congestionDf, aes(x=rd.rp.hour,y=avRep)) + geom_bar(stat="identity")
congestionDf <- congestionDf %>% arrange(desc(avRep))
congestionDf %>% head %>% kable
```

It is shown that 18.00 is the highest congestion hour, and the highest congestion period is from 15.00 to 19.00 (inclusively). For the morning period, a congestion peak at hour 8.00 is detected.

##Congestion days
We now present another distribution to indicate the congestion across the weekdays. A similar approach was followed by filtering out the special reports and questions.
```{r}
congestionDf <- df[df$rd.rp.stid < 6,] %>% group_by(rd.rp.wday) %>% summarize(avRep=mean(rd.rp.stid), numRep=n())
ggplot(congestionDf, aes(x=rd.rp.wday,y=avRep)) + geom_bar(stat="identity")
```

It is observed that the day having the highest overall congestion is Monday. The lowest congestion is on Friday, and the lowest congestion on a working day is on Saturday.

## Congestion regions
Now, with the same filtering technique on report types, we take a look at roads with the highest congestion. A new filtering approach is added here based on the number of reports. This aproach was not needed before since each hour or weekday typically had a significant number of reports. This is not the case here, as some unpopular roads can have a trivial number of reports throughout the whole dataset.
```{r}
congestionDf <- df[df$rd.rp.stid < 6,] %>% group_by(rd.nm) %>% summarize(avRep=mean(rd.rp.stid), numRep=n()) %>% arrange(desc(avRep))
head(congestionDf,1) %>% kable
```

As shown above, the road proposed to have the highest congestion only relies on 3 reports, which is too trivial to consider. To remedy this, we threshold our results according to the mean number of reports among all roads.
```{r}
meanNmRep <- mean(congestionDf$numRep)
congestionDfHighNmRep <- congestionDf[congestionDf$numRep > meanNmRep,]
congestionDfHighNmRep %>% head %>% kable
```

This makes much more sense now.

##Highest user activity hours
Now we examine the hours which have the highest overall user activity, in terms of reports submitted, regardless of the report type.
```{r}
ggplot(df, aes(factor(rd.rp.hour))) + geom_bar()
```

As expected, the two time periods with highest activity are the morning from 7.00 to 10.00 and the period from 16.00 to 19.00, with local maximums at 8.00 & 17.00 respectively.


# Road directions
The `rd.nm` column typically consists of 2 parts: the name of the **road** itself and its **direction**. So far, the studies that were distributed across the roads were with respect to each direction. By dividing this column into two columns, we can now distribute with respect to the road as a whole.
```{r warning=FALSE}
genRoadDf <-  df %>% separate(rd.nm,c("rd.name","rd.dir"),";")
```
Now we shall proceed by getting the **roads** having the highest user activity (& neglecting the uni-directional ones). We will then do some analysis on the directions of a highly active road. Such analysis can lead to observations about this specific road and properties of each direction in it, that were not visible before.
```{r}
genRoadDf[!is.na(genRoadDf$rd.dir),] %>% group_by(rd.name) %>% summarize(numRep=n(),numDir=length(unique(rd.dir))) %>% arrange(desc(numRep)) %>% head %>% kable
```

While "Da2ery" road is the one with the highest activity, it is distributed among 24 directions. So we choose to do our analysis on "Kobry 6 October" road which has only 4 directions, to produce results that are visually clearer.
```{r}
genRoadDf <- genRoadDf[genRoadDf$rd.name=="Kobry 6 October",]
ggplot(genRoadDf, aes(rd.rp.hour)) + geom_bar() + facet_wrap(~ rd.dir)
```

An interesting observation here is that roads whose directions are **to Ta7rir** exhibit higher user activity than roads whose directions are **from Ta7rir**.

# Summary
We have arrived at the end of our analysis on the traffic report data crawled from bey2ollak.com website. We began by cleaning the dataset by getting rid of unwanted rows and columns, and converting the desired ones into more usable forms. Understanding the column meanings was done by querying the dataset and comparing it with bey2ollak.com response data. Interesting insights have been obtained regarding user activity and traffic congestion times. Furhtermore, some observations regarding highly congested roads were made, which agree with our daily struggles on such roads. Lastly, we singled out "Kobry 6 October" road and focused on the user activity on its different directions.