GPS Data Analysis

# Notebook Setup
# import all modules
import pandas as pd
import numpy as np
import seaborn as sns
import matplotlib.pyplot as plt
import geopy.distance
# Ignore warnings
import warnings
warnings.filterwarnings("ignore")
% matplotlib inline

gps_data=pd.read_csv("gps_data.csv", parse_dates=["timestamp"], index_col="timestamp")

# Print the first five rows
gps_data.head()

# Print the last 5 rows
gps_data.tail()

gps_data.info()

gps_data.describe()

sns.set_style("darkgrid")
fig, ax=plt.subplots(figsize=(8,6))
sns.distplot(gps_data["accuracy"])

ax.set_xlabel("Accuracy", fontsize=13)
ax.set_ylabel("differential (accuracy)", fontsize=13)
ax.set_title("Distribution of Accuracy");
ax.set_xlim(left=0, right=500);

# Remove the rows with accuracy level > 57 m
# The values out of the third quartile
gps_data= gps_data[gps_data["accuracy"]< 57]

# See the size and other parameters of the dataset after removal
gps_data.info()
gps_data.describe()

# Add a column showing the time differences between two successive data points in seconds
gps_data["time_difference_sec"]=gps_data["utc_time_sec"].diff()

# Print the details of time_difference column with .describe() method 
print(gps_data.time_difference_sec.describe())

# Print the head 
gps_data.head(3)

fig, ax=plt.subplots(figsize=(8,6))
gps_data.time_difference_sec.plot(marker="o", linewidth=0)

ax.set_xlabel("Week", fontsize=13)
ax.set_ylabel("Time difference in seconds", fontsize=13)
ax.set_title("Time difference between two subsequent data sample");

fig, ax=plt.subplots(figsize=(8,4))
sns.distplot(gps_data["time_difference_sec"].dropna(), 
                 rug=True, 
                 rug_kws={"color":"black"})

ax.set_xlabel("Time difference in seconds", fontsize=13)
ax.set_ylabel("d (time difference)", fontsize=13);
ax.set_title("Distribution of Time differences");
#ax.set_xlim(0, 4000);

# Display the joint distribution of longitudes and lattitudes
g= sns.JointGrid(x="longitude", 
                 y= "latitude", 
                 data=gps_data, 
                 height=7)

g= g.plot_joint(sns.kdeplot, cmap="Blues_d")
g= g.plot(sns.regplot, sns.distplot)

# See how the distance function works

# Two coordinate samples with 0.1 change in lattitudes
coords1 = (20.0, 70)
coords2 = (20.1, 70)

# Pass the coordinates to the distance function and 
# use the `.km` attribute to get the distance in kilometer
print("When the other component is constant")
print("0.1 change in latitude corresponds", geopy.distance.distance(coords1, coords2).km, " km change")

# Another coordinate samples with 0.1 change in longitude
coords2 = (20, 70)
coords3 = (20, 70.1)

print( "0.1 change in longitude corresponds", geopy.distance.distance(coords2, coords3).km, "km change")

# Represent a location as a latitude and longitude tuple and
# create the locations list whose elements are the tuples consist of latitude and longitude values
locations=list(zip(gps_data.latitude, gps_data.longitude))

# Make a list of succesive distances in kms
# by using list comprehensive method and geopy.distance.distance function
distance_lst=[geopy.distance.distance(locations[i-1], locations[i]).km for i in range(1, len(locations))]

# Print the first 3 elements of the distance_lst
distance_lst[:3]

# Since the number of intervals between n succesive points is n-1
# We need to insert a 0 (zero) in the first position of the list in order to 
# match the index size of the gps_data dataframe
distance_lst.insert(0, 0)

# Convert the distance list to an array 
distance_ary=np.array(distance_lst)

# Round the distance values to 4 decimals
distance_ary=np.round(distance_ary, 4)

# Add a new column of distances in kms: distance_km
gps_data["distance_km"]= distance_ary

# See the head of the new data frame
gps_data.head(3)

# Excrat the weekdays from the gps_data: weekdays
weekdays=gps_data[gps_data.index.weekday < 5]

# Excrat the weekends from the gps_data: weekend
weekend=gps_data[gps_data.index.weekday >= 5]

# Define the barplot function: bar_plot
def bar_plot(df, label_x, label_y, str_title):
    '''Takes a dataframe and three strings: label of x-axis, label of y-axis and title
        returns a bar plot'''
    plt.style.use('fivethirtyeight')
    fig, ax=plt.subplots()
    ax= df.plot(kind= "bar", 
                figsize=(10,5), 
                rot=70,
                fontsize=11)        
    ax.set_xlabel(label_x)
    ax.set_ylabel(label_y)
    ax.set_title(str_title)
    plt.show()

# Count the daily data records in weekdays: daily_count
# Groupby the gps_data daily by extracting the "day" from the datetime index and 
# count any of the column (without NaN values)
daily_weekdays_count= gps_data.groupby(gps_data.index.dayofweek).accuracy.count()

# Set the index with the names of days
daily_weekdays_count.index=['Mon', 'Tue', 'Wed', 'Thu', 'Fri', "Sat", "Sun"]

labelx="Day"
labely= "Data counts"
title="Daily data count"
bar_plot(daily_weekdays_count, labelx, labely, title)

# Count the hourly data records in weekdays: hourly_weekdays_count
# Groupby the weekdays hourly by extracting the "hour" from the datetime index and 
# count any of the column (without NaN values)
hourly_weekdays_count= weekdays.groupby(weekdays.index.hour).accuracy.count()

# Set the index with the keys of the group object, here hours
hourly_weekdays_count.index= hourly_weekdays_count.keys()

labelx="Hour"
labely= "Data counts"
title="Weekdays hourly data count"
bar_plot(hourly_weekdays_count, labelx, labely, title)

# Count the hourly data records in weekends: hourly_weekends_count
# Groupby the weekends hourly by extracting the "hour" from the datetime index and 
# count any of the column (without NaN values)
hourly_weekends_count= weekend.groupby(weekend.index.hour).accuracy.count()

# Set the index with the keys of the group object, here hours
hourly_weekends_count.index= hourly_weekends_count.keys()

labelx="Hour"
labely= "Data counts"
title="Weekends hourly data count"
bar_plot(hourly_weekdays_count, labelx, labely, title)

# Slice the period 0-5 hours from the weekdays df: night_weekdays
night_weekdays=weekdays[(weekdays.index.hour>=0) & (weekdays.index.hour<5)]
                        
# Days in which we get data during 0-5 hours
# Extract the unique days by .unique() method 
days0_5=night_weekdays.index.day.unique()

# Print all the weekdays first
print('All the weekdays of the dataset:', "\n", f'{set(weekdays.index.day)}')
print(f'({len(set(weekdays.index.day))} days)', "\n")

# Print only the days we record data
print("Data recording days during 0h-5h:", "\n", list(days0_5))

# Slice the period 5-19 hours from the weekdays df: night_weekdays
daytime_weekdays=weekdays[(weekdays.index.hour>=5) & (weekdays.index.hour<19)]

# Days in which we get data during 5-19 hours
days5_19= daytime_weekdays.index.day.unique()

# Print only the days we record data
print("Data recording days during 5h-19h:", "\n", list(days5_19))
print(f'({len(list(days5_19))} days)')

# Slice the period 19-24 hours from the weekdays df: night_weekdays
evening_weekdays=weekdays[(weekdays.index.hour>=19) & (weekdays.index.hour<24)]

# Days in which we get data during 19_24 hours
days19_24= evening_weekdays.index.day.unique()

# Print only the days data is recorded
print("Data recording days during 19_24:", "\n", list(days19_24))
print(f'({len(list(days19_24))} days)')

# Find the average daily distance_km values in weekdays df: daily_movement
# Groupby the the distance_km column values daily and aggregate by mean
daily_movement=weekdays.groupby(weekdays.index.dayofweek).distance_km.mean()

# Set the index by the names of the weekdays
daily_movement.index = ['Mon', 'Tue', 'Wed', 'Thu', 'Fri']

# Plot the chart
labelx="Day"
labely= "Total movement (km)"
title="Weekdays daily movement"
bar_plot(daily_movement, labelx, labely, title)

def group_plot_coord(groupby_obj, title_func):
    '''Takes a groupby object and a function (to create the titles of the subplots)
    returns a subplot for each of the group'''
    
    # Use the 'fivethirtyeight' style
    plt.style.context('fivethirtyeight')
    # associate the row number argument in plt.subplots with groupby obj
    row_number=int(len(groupby_obj)/2)+1
    
    # Create figure and axes
    # We fix the number of column of the figure by ncols=2 and
    # Let the row_number change by the size of the groupby object
    # We will delete the last empty subplots at the end
    fig, axes = plt.subplots(nrows=row_number, ncols=2,
                             figsize=(16, 16),
                             sharex=True, 
                             sharey=True)

    # Flatten the axes and convert into a list in order to be able to pop the elements one by one
    axes_lst=list(axes.flatten())

    # Iterate over the groups in the groupby object
    for group_name, group in groupby_obj:
        # Assign the first ax in the ax_lst as current ax
        ax=axes_lst.pop(0)
        # Create the scatter plot of latitudes and longitudes of entire gps_data df
        ax.scatter(gps_data["longitude"], gps_data["latitude"],
                   color=".5",
                   edgecolor="white",
                   s=50)
        # Plot the lattidues and longitudes of each group on the same ax
        group.plot(x="longitude", y="latitude", 
                   ax=ax, 
                   kind="scatter", 
                   legend=True,
                   label=group_name, 
                   s=60, 
                   c="r",
                  edgecolor="white")
        # Set the title of each subplot
        ax.set_title(title_func(group_name))
        plt.show()
    
    # Use matplotlib .remove() method to delete the last empty subplot
    [ax.remove() for ax in axes_lst];

# Groupby the weekdays dataframe based on days of the week: daily_group
daily_groups=weekdays.groupby(weekdays.index.dayofweek)

# Create a dictionary of day names in order to use for the titles of subplots
day_dict={0:"Monday", 1:"Tuesday", 2:"Wednesday", 3:"Thursday", 4:"Friday"}

# define the title func:title_func
def daily_title(name):
    return day_dict[name]
# Plot the coordinates of all the weekdays
group_plot_coord(daily_groups, daily_title)

# Display the joint distribution of longitudes and lattitudes
def full_joint_plot(df, title_str):
    '''Takes a dataframe and title string and returns a joint distribution of the coordinates'''
    sns.axes_style("darkgrid")
    g= sns.JointGrid(x="longitude", 
                     y= "latitude", 
                     data=df, 
                     height=7)
    g.plot_joint(sns.kdeplot, cmap="Blues_d")
    sns.kdeplot(df["longitude"], ax=g.ax_marg_x, legend=False)
    sns.kdeplot(df["latitude"], ax=g.ax_marg_y, vertical=True, legend=False)
    g.ax_joint.plot(df["longitude"], df["latitude"], "o", ms=5)
    g.fig.suptitle(title_str, x=0.4, y=0.81);

# Slice the rows of the day monday-to-thursday from gps_data: mon_to_thu
mon_to_thu= gps_data[gps_data.index.dayofweek < 4]
title_str="Coordinates distribution monday-to-thursday"
full_joint_plot(mon_to_thu, title_str)

# Slice the rows on fridays from gps_data: fridays
fridays=gps_data[gps_data.index.dayofweek==4]

title_str="Coordinate distribution on fridays"
full_joint_plot(fridays, title_str)

# Find the average hourly distance_km values in mon_to_thu df: hourly_movement1_4
# Groupby the the distance_km column values hourly and aggregate by mean
hourly_movement1_4=mon_to_thu.groupby(mon_to_thu.index.hour).distance_km.mean()

# Set the index by the keys of the group object, here hours
hourly_movement1_4.index = hourly_movement1_4.keys()

# Plot the chart
labelx="Hour"
labely= "Average movement (km)"
title="Monday-to-thursday hourly movement"
bar_plot(hourly_movement1_4, labelx, labely, title)

# drop the outlier
fridays=fridays.drop(pd.to_datetime("2015-03-20 07:21:58"))

# Find the average hourly distance_km values in fridays df: hourly_movement4
# Groupby the the distance_km column values hourly and aggregate by mean
hourly_movement4=fridays.groupby(fridays.index.hour).distance_km.mean()

# Set the index by the keys of the group object, here hours
hourly_movement4.index = hourly_movement4.keys()

# Plot the chart
labelx="Hour"
labely= "Average movement (km)"
title="Fridays hourly movement"
bar_plot(hourly_movement4, labelx, labely, title)

# Slice the morning times from mon_to_thu df
morning0_3=mon_to_thu[(mon_to_thu.index.hour>=5)& (mon_to_thu.index.hour<12)]

# Groupby the morning_weekdays hourly: hourly_morning
hourly_morning0_3=morning0_3.groupby(morning0_3.index.hour)

def hourly_title(name):
    return f'Coordinates during {name}h-{name+1}h'

group_plot_coord(hourly_morning0_3, hourly_title)

# Slice the morning times from mon_to_thu df
afternoon0_3=mon_to_thu[(mon_to_thu.index.hour>=15)& (mon_to_thu.index.hour<19)]

# Groupby the morning_weekdays hourly: hourly_morning
hourly_afternoon0_3=afternoon0_3.groupby(afternoon0_3.index.hour)

def hourly_title(name):
    return f'Coordinates during {name}h-{name+1}h'

group_plot_coord(hourly_afternoon0_3, hourly_title)

# Slice the morning times from fridays df
morning_fridays=fridays[(fridays.index.hour>=6)& (fridays.index.hour<12)]

# Groupby the morning_fridays hourly: hourly_morning_friday
hourly_morning_friday=morning_fridays.groupby(morning_fridays.index.hour)

def hourly_title(name):
    return f'Coordinates during {name}h-{name+1}h'

group_plot_coord(hourly_morning_friday, hourly_title)

# Slice the morning times from fridays df
afternoon_fridays=fridays[(fridays.index.hour>=15)& (fridays.index.hour<24)]

# Groupby the morning_fridays hourly: hourly_morning_friday
hourly_afternoon_friday=afternoon_fridays.groupby(afternoon_fridays.index.hour)

def hourly_title(name):
    return f'Coordinates during {name}h-{name+1}h'

group_plot_coord(hourly_afternoon_friday, hourly_title)

title= "Weekdays morning distribution"
full_joint_plot(morning0_3, title)

title= "Weekdays morning distribution"
full_joint_plot(morning_fridays, title)

# Subset the gps_data df by just taking distance_km column: heatmap_df
# Since we will do summation and averaging we can replace the NaN values with 0 (zero)
heatmap_df=mon_to_thu[["distance_km"]]

# Add "day" column by extracting the day from the datetime index
heatmap_df["day"]= heatmap_df.index.dayofweek

# Add "hour" column by extracting the hour from the datetime index
heatmap_df["hour"]= heatmap_df.index.hour

# Drop the timestamp index
heatmap_df=heatmap_df.reset_index(drop=True)

# Because of the dublicated values we cannot pivot the dataframe directly for the heatmap
# First we need to create a pivot table
p_table = pd.pivot_table(heatmap_df, 
                       values='distance_km', 
                       index=['day', 'hour'], 
                       aggfunc=np.mean)

# Reset the indexes as columns
p_table=p_table.reset_index(level='hour').reset_index()

# Now we can create the pivot dataframe
heatmap_df= p_table.pivot(index="hour", 
                          columns="day", 
                          values="distance_km")

# Display the heatmap
fig, ax = plt.subplots(figsize=(10, 8)) 
ax = sns.heatmap(heatmap_df,
                 cmap= 'Oranges', 
                 xticklabels=["Mon", "Tue", "Wed", "Thu"])

ax.set_xlabel("Day")
ax.set_ylabel("Hour")
ax.set_title("Movement based on hours of from Mon-to-Thu");