Repository: mie-lab/trackintel
Pull request #103 (2021-02-10 11:04:02)

Feat: simple coarse mode prediction
adapted from #65 
added cases for different CRS:
- WGS84 -> haversine distance (with UserWarning)
- None -> WGS84 is assumed (with UserWarning)
- geographic coordinate system (other than WGS84) -> raises UserWarning, no calculation done
- other (projected) system -> geoSeries.length is used

modified test accordingly

Is the setting of the warnings etc. ok like this? 


Patch:
diff --git a/CONTRIBUTING.md b/CONTRIBUTING.md
new file mode 100644
index 0000000..ebf9f05
--- /dev/null
+++ b/CONTRIBUTING.md
@@ -0,0 +1,43 @@
+# How to contribute
+This is a place to collect conventions we agreed upon until we find the right place in the doc for them
+
+## Coding conventions
+### Time stamps
+All timestamps are timezone aware pandas `datetime64[ns, UTC]` objects. The default timezone should be `UTC` but the user should be free to use a different one if he wants. See [Issue 101](https://github.com/mie-lab/trackintel/issues/18). 
+
+### Organization of tests
+See [Issue 23](https://github.com/mie-lab/trackintel/issues/23)
+- The test folder copies the folder structure that the trackintel.trackintel folder has.
+- Every python module has a single test file
+- Every function has 1 test class
+- Every method of this function should test a single property
+
+### Warnings
+There should be no warnings when running tests.
+
+### Integrety of input data
+Functions should never change the input dataframe but rather return an altered copy.
+
+### Adressing geometry columns
+Geometry columns should never be adressed by name but by the geometry attribute [Issue 14](https://github.com/mie-lab/trackintel/issues/14), [Issue 15](https://github.com/mie-lab/trackintel/issues/15)
+
+### Order of functions in code
+The main function should be on the top of the file, the internal/secondary functions should be at the end of the file
+
+### ID management
+All trackintel objects have an ID that is the index of the dataframe [Issue 97](https://github.com/mie-lab/trackintel/issues/97)
+
+### Others
+- We limit all lines to a maximum of 120 characters.
+
+## Version release checklist
+Before you release a new version you should check/modify the following files:
+
+- trackintel/__ version__.py -> update the version number
+
+- docs/index.rst -> update the version number
+
+- check setup.py if new dependencies need greater python version
+
+After the release is tagged, __no new release with the same version number__ can be published on pypi!
+
diff --git a/Readme_developers.md b/Readme_developers.md
deleted file mode 100644
index 68342d1..0000000
--- a/Readme_developers.md
+++ /dev/null
@@ -1,14 +0,0 @@
-## trackintel: version release checklist
-
-Before you release a new version you should check/modify the following files:
-
-
-
-- trackintel/__ version__.py -> update the version number
-
-- docs/index.rst -> update the version number
-
-- check setup.py if new dependencies need greater python version
-
-After the release is tagged, __no new release with the same version number__ can be published on pypi!
-
diff --git a/docs/modules/analysis.rst b/docs/modules/analysis.rst
index 5d327a0..5b4afb3 100644
--- a/docs/modules/analysis.rst
+++ b/docs/modules/analysis.rst
@@ -1,4 +1,10 @@
 Analysis
 ********
 
-The analysis module documentation is yet to be written :).
\ No newline at end of file
+The analysis module provides a variety of functions to analyze the mobility
+datasets computed by trackintel. 
+
+Triplegs
+========
+
+.. autofunction:: trackintel.analysis.transport_mode_identification.predict_transport_mode
diff --git a/trackintel/__init__.py b/trackintel/__init__.py
index 6d5953e..fbc63a7 100644
--- a/trackintel/__init__.py
+++ b/trackintel/__init__.py
@@ -1,3 +1,4 @@
 from .core import *
 from .model import *
-from .io import *
\ No newline at end of file
+from .io import *
+from .analysis import *
\ No newline at end of file
diff --git a/trackintel/analysis/__init__.py b/trackintel/analysis/__init__.py
index e69de29..849b37f 100644
--- a/trackintel/analysis/__init__.py
+++ b/trackintel/analysis/__init__.py
@@ -0,0 +1,1 @@
+from .transport_mode_identification import predict_transport_mode
diff --git a/trackintel/analysis/transport_mode_identification.py b/trackintel/analysis/transport_mode_identification.py
new file mode 100644
index 0000000..7015f80
--- /dev/null
+++ b/trackintel/analysis/transport_mode_identification.py
@@ -0,0 +1,136 @@
+
+import trackintel as ti
+from trackintel.geogr.distances import haversine_dist
+import geopandas as gpd
+import warnings
+import numpy as np
+
+
+def predict_transport_mode(triplegs, method='simple-coarse', **kwargs):
+    """
+    Predicts/imputes the transport mode that was likely chosen to cover the \
+    given tripleg, e.g., car, bicycle, or walk.
+
+    Parameters
+    ----------
+    method: str, {'simple-coarse'}, default 'simple-coarse'
+        The following methods are available for transport mode inference/prediction:
+        'simple-coarse' : Uses simple heuristics to predict coarse transport classes.
+            These include ``{'slow_mobility', 'motorized_mobility', 'fast_mobility'}``.
+            In the default classification, ``slow_mobility`` (<15 km/h) includes transport modes such as
+            walking or cycling, ``motorized_mobility`` (<100 km/h) modes such as car or train, and
+            ``fast_mobility`` (>100 km/h) modes such as high-speed rail or airplanes.
+            These categories are default values and can be overwritten using the keyword argument categories.
+    """
+    if method == 'simple-coarse':
+        # implemented as keyword argument if later other methods that don't use categories are added
+        categories = kwargs.pop('categories', {15/3.6: 'slow_mobility', 100/3.6: 'motorized_mobility',
+                                               np.inf: 'fast_mobility'})
+
+        return predict_transport_mode_simple_coarse(triplegs, categories)
+    else:
+        raise NameError(f'Method {method} not known for predicting tripleg transport modes.')
+
+
+def predict_transport_mode_simple_coarse(triplegs, categories):
+    """
+    Predicts a transport mode out of three coarse classes. Implements a simple speed based heuristic
+     (over the whole tripleg). As such, it is very fast, but also very simple and coarse.
+
+    Parameters
+    ----------
+    triplegs : trackintel triplegs GeoDataFrame
+        The triplegs for the transport mode prediction.
+    categories : dict, optional
+        The categories for the speed classification {upper_boundary:'category_name'}. The unit for the upper boundary
+        is m/s.
+        The default is {15/3.6: 'slow_mobility', 100/3.6: 'motorized_mobility', np.inf: 'fast_mobility'}.
+
+    Raises
+    ------
+    ValueError
+        In case the boundaries of the categories are not in ascending order.
+
+    Returns
+    -------
+    triplegs : trackintel triplegs GeoDataFrame
+        the triplegs with added column mode, containing the predicted transport modes.
+
+    For additional documentation, see
+    :func:`trackintel.analysis.transport_mode_identification.predict_transport_mode`.
+
+    """
+    if not(check_categories(categories)):
+        raise ValueError('the catecories must be in increasing order')
+
+    triplegs = triplegs.copy()
+    wgs = False
+
+    if triplegs.crs == 4326:
+        wgs = True
+
+    elif triplegs.crs is None:
+        wgs = True
+        warnings.warn('Your data is not projected. WGS84 is assumed and for length calculation the haversine '
+                      'distance is used')
+
+    elif triplegs.crs.is_geographic:
+        raise UserWarning('Your data is in a geographic coordinate system, length calculation fails')
+
+    def identify_mode(tripleg, wgs, categories):
+        """
+        Identify the mode based on the (overall) tripleg speed.
+
+        Parameters
+        ----------
+        tripleg : trackintel triplegs GeoDataFrame
+            the tripleg to analyse
+        wgs : bool
+            whether the tripleg is in WGS84 or not.
+        categories : dict
+            the upper boundaries (as keys) and the names of the categories as values.
+
+        Returns
+        -------
+        str
+            the identified mode.
+        """
+        # Computes distance over whole tripleg geometry (using the Haversine distance).
+        if wgs:
+            distance = sum([haversine_dist(pt1[0], pt1[1], pt2[0], pt2[1]) for pt1, pt2
+                            in zip(tripleg.geom.coords[:-1], tripleg.geom.coords[1:])])
+        else:
+            distance = tripleg.geom.length
+
+        duration = (tripleg['finished_at'] - tripleg['started_at']).total_seconds()
+        speed = distance / duration  # The unit of the speed is m/s
+
+        for bound in categories:
+            if speed < bound:
+                return categories[bound]
+
+    triplegs['mode'] = triplegs.apply(lambda l: identify_mode(l, wgs, categories), axis=1)
+    return triplegs
+
+
+def check_categories(cat):
+    """
+    Check if the keys of a dictionary are in ascending order.
+
+    Parameters
+    ----------
+    cat : disct
+        the dictionary to be checked.
+
+    Returns
+    -------
+    correct : bool
+        True if dict keys are in ascending order False otherwise.
+
+    """
+    correct = True
+    bounds = list(cat.keys())
+    for i in range(len(bounds)-1):
+        if bounds[i] >= bounds[i+1]:
+            correct = False
+    return correct
diff --git a/trackintel/geogr/distances.py b/trackintel/geogr/distances.py
index 3c6df0c..74848e6 100644
--- a/trackintel/geogr/distances.py
+++ b/trackintel/geogr/distances.py
@@ -8,6 +8,7 @@ from sklearn.metrics import pairwise_distances
 
 from trackintel.geogr.point_distances import haversine_dist
 from trackintel.geogr.trajectory_distances import dtw, frechet_dist
+import warnings
 
 
 def calculate_distance_matrix(X, Y=None, dist_metric='haversine', n_jobs=0, **kwds):
@@ -161,4 +162,4 @@ def meters_to_decimal_degrees(meters, latitude):
     float
         An approximation of a distance (given in meters) in degrees.
     """
-    return meters / (111.32 * 1000.0 * cos(latitude * (pi / 180.0)))
\ No newline at end of file
+    return meters / (111.32 * 1000.0 * cos(latitude * (pi / 180.0)))
diff --git a/trackintel/io/file.py b/trackintel/io/file.py
index 9e526a4..0d06374 100644
--- a/trackintel/io/file.py
+++ b/trackintel/io/file.py
@@ -166,9 +166,9 @@ def write_triplegs_csv(triplegs, filename, *args, **kwargs):
     filename : str
         The file to write to.
     """
-
-    gdf = triplegs.copy()
-    gdf[gdf.geometry.name] = triplegs.geometry.apply(wkt.dumps)
+    geo_col_name = triplegs.geometry.name
+    gdf = pd.DataFrame(triplegs, copy=True)
+    gdf[geo_col_name] = triplegs.geometry.apply(wkt.dumps)
     gdf.to_csv(filename, index=True, *args, **kwargs)
 
 
@@ -231,8 +231,9 @@ def write_staypoints_csv(staypoints, filename, *args, **kwargs):
     filename : str
         The file to write to.
     """
-    gdf = staypoints.copy()
-    gdf[gdf.geometry.name] = gdf.geometry.apply(wkt.dumps)
+    geo_col_name = staypoints.geometry.name
+    gdf = pd.DataFrame(staypoints, copy=True)
+    gdf[geo_col_name] = staypoints.geometry.apply(wkt.dumps)
     gdf.to_csv(filename, index=True, *args, **kwargs)
 
 
@@ -281,7 +282,7 @@ def write_locations_csv(locations, filename, *args, **kwargs):
     filename : str
         The file to write to.
     """
-    gdf = locations.copy()
+    gdf = pd.DataFrame(locations, copy=True)
     gdf['center'] = locations['center'].apply(wkt.dumps)
     if 'extent' in gdf.columns:
         gdf['extent'] = locations['extent'].apply(wkt.dumps)
diff --git a/trackintel/model/triplegs.py b/trackintel/model/triplegs.py
index 69f1206..ee077dd 100644
--- a/trackintel/model/triplegs.py
+++ b/trackintel/model/triplegs.py
@@ -82,3 +82,9 @@ class TriplegsAccessor(object):
         """Filter triplegs with a geo extent.
         See :func:`trackintel.preprocessing.filter.spatial_filter`."""
         return ti.preprocessing.filter.spatial_filter(self._obj, *args, **kwargs)
+
+    def predict_transport_mode(self, *args, **kwargs):
+        """Predict/impute the transport mode with which each tripleg was likely covered.
+        See :func:`trackintel.analysis.transport_mode_identification.predict_transport_mode`.
+        """
+        return ti.analysis.transport_mode_identification.predict_transport_mode(self._obj, *args, **kwargs)
diff --git a/trackintel/preprocessing/staypoints.py b/trackintel/preprocessing/staypoints.py
index 70d8e10..1846f33 100644
--- a/trackintel/preprocessing/staypoints.py
+++ b/trackintel/preprocessing/staypoints.py
@@ -56,10 +56,11 @@ def generate_locations(staypoints,
     
     if agg_level not in ['user', 'dataset']:
         raise AttributeError("The parameter agg_level must be one of ['user', 'dataset'].")
+    if method not in ['dbscan']:
+        raise AttributeError("The parameter method must be one of ['dbscan'].")
     
     # initialize the return GeoDataFrames
     ret_stps = staypoints.copy()
-    ret_loc = gpd.GeoDataFrame([], columns=['user_id', 'id', 'center', 'extent'])  # todo: define default dtypes
     
     if method=='dbscan':
 
@@ -74,6 +75,7 @@ def generate_locations(staypoints,
             
         if agg_level == 'user':
             location_id_counter = 0
+            # TODO: change into groupby
             for user_id_this in ret_stps["user_id"].unique():
                 # Slice staypoints array by user. This is not a copy!
                 user_staypoints = ret_stps[ret_stps["user_id"] == user_id_this]  
@@ -120,33 +122,38 @@ def generate_locations(staypoints,
         # filter outlier
         ret_loc = ret_loc.loc[ret_loc['location_id'] != -1]
         
+        ret_loc['center'] = None # initialize
         # locations with only one staypoints is of type "Point"
         point_idx = ret_loc.geom_type == 'Point'
-        ret_loc['center'] = 0 # initialize
-        ret_loc.loc[point_idx, 'center'] = ret_loc.loc[point_idx, 'geom']
+        if not ret_loc.loc[point_idx].empty:
+            ret_loc.loc[point_idx, 'center'] = ret_loc.loc[point_idx, 'geom']
         # locations with multiple staypoints is of type "MultiPoint"
-        ret_loc.loc[~point_idx, 'center'] = ret_loc.loc[~point_idx, 'geom'].apply(
-            lambda p: Point(np.array(p)[:,0].mean(), np.array(p)[:,1].mean()))
+        if not ret_loc.loc[~point_idx].empty:
+            ret_loc.loc[~point_idx, 'center'] = ret_loc.loc[~point_idx, 'geom'].apply(
+                lambda p: Point(np.array(p)[:,0].mean(), np.array(p)[:,1].mean()))
         
         # extent is the convex hull of the geometry
-        ret_loc['extent'] = ret_loc['geom'].apply(lambda p: p.convex_hull)                                                                                        
-        # convex_hull of one point would be a Point and two points a Linestring, 
-        # we change them into Polygon by creating a buffer of epsilon around them.
-        pointLine_idx = (ret_loc['extent'].geom_type == 'LineString') | (ret_loc['extent'].geom_type == 'Point')
+        ret_loc['extent'] = None # initialize
+        if not ret_loc.empty:
+            ret_loc['extent'] = ret_loc['geom'].apply(lambda p: p.convex_hull)                                                                                        
+            # convex_hull of one point would be a Point and two points a Linestring, 
+            # we change them into Polygon by creating a buffer of epsilon around them.
+            pointLine_idx = (ret_loc['extent'].geom_type == 'LineString') | (ret_loc['extent'].geom_type == 'Point')
         
-        # Perform meter to decimal conversion if the distance metric is haversine
-        if distance_matrix_metric == 'haversine':
-            ret_loc.loc[pointLine_idx, 'extent'] = ret_loc.loc[pointLine_idx].apply(
-                lambda p: p['extent'].buffer(meters_to_decimal_degrees(epsilon, p['center'].y)), axis=1)
-        else:
-            ret_loc.loc[pointLine_idx, 'extent'] = ret_loc.loc[pointLine_idx].apply(
-                lambda p: p['extent'].buffer(epsilon), axis=1)
+            if not ret_loc.loc[pointLine_idx].empty:
+                # Perform meter to decimal conversion if the distance metric is haversine
+                if distance_matrix_metric == 'haversine':
+                    ret_loc.loc[pointLine_idx, 'extent'] = ret_loc.loc[pointLine_idx].apply(
+                        lambda p: p['extent'].buffer(meters_to_decimal_degrees(epsilon, p['center'].y)), axis=1)
+                else:
+                    ret_loc.loc[pointLine_idx, 'extent'] = ret_loc.loc[pointLine_idx].apply(
+                        lambda p: p['extent'].buffer(epsilon), axis=1)
         
         ret_loc = ret_loc.set_geometry('center')
         ret_loc = ret_loc[['user_id', 'location_id', 'center', 'extent']]
         ret_loc.rename(columns={'location_id':'id'}, inplace=True)
-        
-    ## ensure dtype consistency 
+
+    # ensure dtype consistency
     # locs id (generated by this function) should be int64
     ret_loc['id'] = ret_loc['id'].astype('int64')
     ret_stps['location_id'] = ret_stps['location_id'].astype('int64')
@@ -156,7 +163,7 @@ def generate_locations(staypoints,
     return ret_stps, ret_loc
 
 
-def create_activity_flag(staypoints, method='time_threshold', time_threshold=5, activity_column_name='activity'):
+def create_activity_flag(staypoints, method='time_threshold', time_threshold=5, activity_column_name='activity'):           
     """
 
     Parameters


Test patch:
diff --git a/tests/analysis/test_transport_mode_identification.py b/tests/analysis/test_transport_mode_identification.py
new file mode 100644
index 0000000..2e25f65
--- /dev/null
+++ b/tests/analysis/test_transport_mode_identification.py
@@ -0,0 +1,74 @@
+# -*- coding: utf-8 -*-
+
+import datetime
+import os
+
+import pytest
+import numpy as np
+import geopandas as gpd
+from numpy.testing import assert_almost_equal
+
+import trackintel as ti
+from trackintel.io.dataset_reader import read_geolife
+
+
+class TestTransportModeIdentification:
+    def test_check_empty_dataframe(self):
+        """Assert that the method does not work for empty DataFrames (but
+        that the rest works fine, e.g., method signature)."""
+        tpls = ti.read_triplegs_csv('tests/data/triplegs_transport_mode_identification.csv', sep=';')
+        empty_frame = tpls[0:0]
+        with pytest.raises(AssertionError):
+            empty_frame.as_triplegs.predict_transport_mode(method='simple-coarse')
+
+    def test_simple_coarse_identification_no_crs(self):
+        """Assert that the simple-coarse transport mode identification throws the correct \
+warning and and yields the correct results for WGS84 ."""
+        tpls = ti.read_triplegs_csv('tests/data/triplegs_transport_mode_identification.csv', sep=';')
+
+        with pytest.warns(UserWarning, match='Your data is not projected. WGS84 is assumed and for length \
+calculation the haversine distance is used'):
+            tpls_transport_mode = tpls.as_triplegs.predict_transport_mode(method='simple-coarse')
+
+        assert tpls_transport_mode.iloc[0]['mode'] == 'slow_mobility'
+        assert tpls_transport_mode.iloc[1]['mode'] == 'motorized_mobility'
+        assert tpls_transport_mode.iloc[2]['mode'] == 'fast_mobility'
+
+    def test_simple_coarse_identification_wgs_84(self):
+        """Asserts the correct behaviour with data in wgs84"""
+        tpls = ti.read_triplegs_csv('tests/data/triplegs_transport_mode_identification.csv', sep=';')
+        tpls_2 = tpls.set_crs(epsg=4326)
+        tpls_transport_mode_2 = tpls_2.as_triplegs.predict_transport_mode(method='simple-coarse')
+
+        assert tpls_transport_mode_2.iloc[0]['mode'] == 'slow_mobility'
+        assert tpls_transport_mode_2.iloc[1]['mode'] == 'motorized_mobility'
+        assert tpls_transport_mode_2.iloc[2]['mode'] == 'fast_mobility'
+
+    def test_simple_coarse_identification_projected(self):
+        """Asserts the correct behaviour with data in projected coordinate systems"""
+        tpls = ti.read_triplegs_csv('tests/data/triplegs_transport_mode_identification.csv', sep=';')
+        tpls_2 = tpls.set_crs(epsg=4326)
+        tpls_3 = tpls_2.to_crs(epsg=2056)
+        tpls_transport_mode_3 = tpls_3.as_triplegs.predict_transport_mode(method='simple-coarse')
+        assert tpls_transport_mode_3.iloc[0]['mode'] == 'slow_mobility'
+        assert tpls_transport_mode_3.iloc[1]['mode'] == 'motorized_mobility'
+        assert tpls_transport_mode_3.iloc[2]['mode'] == 'fast_mobility'
+
+    def test_simple_coarse_identification_geographic(self):
+        """asserts the correct behaviour with data in geographic coordinate systems"""
+        tpls = ti.read_triplegs_csv('tests/data/triplegs_transport_mode_identification.csv', sep=';')
+        tpls_2 = tpls.set_crs(epsg=4326)
+        tpls_4 = tpls_2.to_crs(epsg=4269)
+        with pytest.raises(UserWarning,
+                           match='Your data is in a geographic coordinate system, length calculation fails'):
+            tpls_4.as_triplegs.predict_transport_mode(method='simple-coarse')
+
+    def test_check_categories(self):
+        """asserts the correct identification of valid category dictionaries"""
+        tpls = ti.read_triplegs_csv('tests/data/triplegs_transport_mode_identification.csv', sep=';')
+        correct_dict = {2: 'cat1', 7: 'cat2', np.inf: 'cat3'}
+
+        assert ti.analysis.transport_mode_identification.check_categories(correct_dict)
+        with pytest.raises(ValueError):
+            incorrect_dict = {10: 'cat1', 5: 'cat2', np.inf: 'cat3'}
+            tpls.as_triplegs.predict_transport_mode(method='simple-coarse', categories=incorrect_dict)
diff --git a/tests/data/triplegs_transport_mode_identification.csv b/tests/data/triplegs_transport_mode_identification.csv
new file mode 100644
index 0000000..3006671
--- /dev/null
+++ b/tests/data/triplegs_transport_mode_identification.csv
@@ -0,0 +1,4 @@
+;user_id;started_at;finished_at;geom
+0;1;2015-11-27T08:00:00Z;2015-11-27T10:00:00Z;LINESTRING (8.5487794879999992 47.3765150500000018, 8.5276999999999994 47.3993499999999983, 8.5116899999999998 47.4085400000000021)
+1;1;2015-11-27T12:00:00Z;2015-11-27T14:00:00Z;LINESTRING (8.5633999999999997 47.9560000000000031, 8.6456000000000000 47.2334499999999977, 8.4534000000000002 47.2633999999999972)
+2;1;2015-11-27T14:00:00Z;2015-11-27T14:30:00Z;LINESTRING (8.5633999999999997 47.9560000000000031, 8.6456000000000000 47.2334499999999977, 8.4534000000000002 47.2633999999999972)