From 1cc52aef8aac5d7f33c38ccf49c235ffc7778f2e Mon Sep 17 00:00:00 2001 From: Chance Haycock Date: Thu, 13 Aug 2020 11:28:44 +0100 Subject: [PATCH 1/5] Extra ScanStatistic Functionality --- SpatialScan/scanstatistic.py | 33 +++++++++++++++++++++++++++------ 1 file changed, 27 insertions(+), 6 deletions(-) diff --git a/SpatialScan/scanstatistic.py b/SpatialScan/scanstatistic.py index 26a7784..83c7364 100644 --- a/SpatialScan/scanstatistic.py +++ b/SpatialScan/scanstatistic.py @@ -104,22 +104,25 @@ def plot(self, metric="av_lhd_score_EBP"): else: logging.info(" Results not populated. Call `run()` first.") - def highest_region(self): + def highest_region(self, metric='l_score_EBP'): """Return highest region""" if isinstance(self.all_results, pd.DataFrame): - return self.all_results.iloc[0] + data = self.all_results.sort_values(metric, ascending=False) + return data.iloc[0] logging.info("Results not populated. Call `run()` first.") - def plot_region_time_series(self, rank=0, legend=False): + def plot_region_time_series(self, rank=0, legend=False, metric='l_score_EBP'): if not isinstance(self.all_results, pd.DataFrame): raise TypeError('Run the scan first') - region = make_region_from_res(self.all_results, rank=rank) + data = self.all_results.sort_values(metric, ascending=False) + region = make_region_from_res(data, rank=rank) plot_region_time_series(region, self.forecast, add_legend=legend) - def plot_region_by_rank(self, rank=0, legend=False): + def plot_region_by_rank(self, rank=0, legend=False, metric='l_score_EBP'): if not isinstance(self.all_results, pd.DataFrame): raise TypeError('Run the scan first') - plot_region_by_rank(rank, self.all_results, self.forecast, add_legend=legend) + data = self.all_results.sort_values(metric, ascending=False) + plot_region_by_rank(rank, data, self.forecast, add_legend=legend) def model_settings(self): settings = self.__dict__.copy() @@ -129,3 +132,21 @@ def model_settings(self): del settings['all_results'] del settings['grid_results'] print(settings) + + def rerun_forecast(self): + self.forecast = count_baseline( + self.processed, + self.days_in_past, + self.days_in_future, + self.ts_method, + alpha=self.alpha, + beta=self.beta, + gamma=self.gamma, + kern=self.kernel, + ) + + def rerun_scan(self): + # Assumes everything remains the same up to scanning + print('Using cached processed and forecast data to rebuild scan') + self.all_results = scan(self.forecast, self.grid_resolution) + self.grid_results = database_results(self.all_results) From cc358f320c0bc2e5ecf8f403ee97e9f4563c0f6b Mon Sep 17 00:00:00 2001 From: Chance Haycock Date: Thu, 13 Aug 2020 11:29:43 +0100 Subject: [PATCH 2/5] Added Network Scan --- SpatialScan/network/__init__.py | 0 SpatialScan/network/network_scan.py | 305 ++++++++++++++++++++++++++ SpatialScan/network/path_generator.py | 112 ++++++++++ 3 files changed, 417 insertions(+) create mode 100644 SpatialScan/network/__init__.py create mode 100644 SpatialScan/network/network_scan.py create mode 100644 SpatialScan/network/path_generator.py diff --git a/SpatialScan/network/__init__.py b/SpatialScan/network/__init__.py new file mode 100644 index 0000000..e69de29 diff --git a/SpatialScan/network/network_scan.py b/SpatialScan/network/network_scan.py new file mode 100644 index 0000000..8552fe9 --- /dev/null +++ b/SpatialScan/network/network_scan.py @@ -0,0 +1,305 @@ +"""Network Scan functionality""" +import logging +import os + +import geopandas as gpd +import numpy as np +import pandas as pd +import osmnx as ox +from shapely.geometry import Point +import networkx as nx + +from SpatialScan.network.path_generator import PathGenerator +from SpatialScan.likelihood import likelihood_ratio + + +class NetworkScan(PathGenerator): + + """Functionality to carry out Scan on Road Networks""" + def __init__( + self, + forecast, + borough="Westminster", + min_path_length=50, + max_path_length=500, + **kwargs + ): + self.forecast = forecast + self.borough = borough + self.borough_file_path = "../../data/ESRI/London_Borough_Excluding_MHW.shp" + + logging.info("Setting up Network Scan") + self.borough_polygon = get_borough_polygon(self.borough, self.borough_file_path) + self.network = get_network_from_polygon(self.borough_polygon) + logging.info("Restricting readings to network") + self.network_forecast = restrict_readings_to_network( + self.network, self.forecast + ) + detector_edges = self.get_detector_edges() + super().__init__( + self.network, min_path_length, max_path_length, detector_edges, **kwargs + ) + logging.info("Calculating paths of lengths between %f and %f", self.min_path_length, self.max_path_length) + self.generate_paths() + + self.t_min = forecast["measurement_start_utc"].min() + self.t_max = forecast["measurement_end_utc"].max() + self.agg_df = None + self.results = None + self.agg_results = None + logging.info("Setup complete.") + + def get_detector_edges(self): + """Return edges of network with active detectors""" + detector_edges = self.network_forecast.drop_duplicates( + subset=["nearest_u", "nearest_v"] + )[["nearest_u", "nearest_v"]] + + detector_edges = detector_edges.astype(int) + + return list(zip(detector_edges.nearest_u, detector_edges.nearest_v)) + + def aggregate_to_edges(self): + """Sum vehicle counts per edge""" + self.network, self.agg_df = aggregate_edge_event_count( + self.network, self.network_forecast, self.detector_edges + ) + + def scan(self): + """Carry out Network Scan, and calculte metrics""" + + t_ticks = pd.date_range(start=self.t_min, end=self.t_max, freq="H") + + return_dict = {} + path_count = 0 + + for time, tick in enumerate(t_ticks[:-1]): + for source in self.network_paths.keys(): + targets = self.network_paths[source].keys() + for target in targets: + + base, count = path_event_count( + self.network, + self.network_paths[source][target], + time, + len(t_ticks) + 1, + ) + + likelihood = likelihood_ratio(base, count) + + return_dict[path_count] = { + "source": source, + "target": target, + "base": base, + "count": count, + "score": likelihood, + "path": self.network_paths[source][target], + "path_graph": self.network.subgraph( + nodes=self.network_paths[source][target] + ), + "measurement_start_utc": tick, + "measurement_end_utc": self.t_max, + } + path_count += 1 + + if time % 8 == 0: + logging.info("Scan progress: %.2f%%", (time + 1) * 100 / len(t_ticks)) + + logging.info("%d space-time paths searched", path_count) + + self.results = pd.DataFrame.from_dict(return_dict, "index") + + def highest_path(self): + """Path with highest EBP score""" + return self.results.iloc[0]["path"] + + + def aggregate_results_to_edges(self): + """Calculate mean EBP score per edge""" + + t_ticks = pd.date_range(start=self.t_min, end=self.t_max, freq="H") + + agg_network = self.network.copy() + + num_aggs = 0 + + # Set Defaults + for time in t_ticks[:-1]: + nx.set_edge_attributes(agg_network, 1, "score_t={}".format(time)) + + agg_dict = {} + for i, (source, target) in enumerate(self.detector_edges): + + spatial_cond = (self.results["path"].astype(str).str.contains("{}, {}".format(source, target)))\ + | (self.results["path"].astype(str).str.contains("{}, {}".format(target, source))) + + for tick in t_ticks[:-1]: + + cond = spatial_cond & (self.results["measurement_start_utc"] <= tick) + sub_df = self.results[cond] + + if sub_df.empty: # TODO - if done properly this shouldnt happen + continue + + score = sub_df["score"].mean() + + agg_dict[num_aggs] = { + "source": source, + "target": target, + "measurement_start_utc": tick, + "measurement_end_utc": self.t_max, + "score": score, + } + num_aggs += 1 + + nx.set_edge_attributes( + agg_network, {(source, target, 0): score}, name="score_t={}".format(tick) + ) + + if i % 50 == 0: + logging.info("Aggregation progress: %.2f%%", i * 100 / len(self.detector_edges)) + + self.agg_results = pd.DataFrame.from_dict(agg_dict, "index") + self.network = agg_network + + +def get_borough_polygon( + borough="Westminster", + boroughs_file="../../data/ESRI/London_Borough_Excluding_MHW.shp", +): + """Fetch borough polygon information""" + + if not os.path.isfile(boroughs_file): + raise ValueError('File does not exist') + + london_df = gpd.read_file(boroughs_file) + london_df = london_df.to_crs(epsg=4326) + return london_df[london_df["NAME"] == borough].iloc[0]["geometry"] + + +def get_network_from_polygon(polygon): + """Fetch network within polygon bounds from osmnx""" + + # Roads of Interest + roi = '["highway"~"motorway|motorway_link|primary|' + roi += 'primary_link|secondary|secondary_link|trunk|trunk_link|tertiary|tertiary_link"]' + network = ox.graph_from_polygon( + polygon, network_type="drive", simplify=True, custom_filter=roi + ) + return nx.MultiGraph(network) # Note - not directional + + +def restrict_readings_to_polygon(forecast, polygon): + + """"Keep readings which lie within polygon""" + + detectors = forecast.drop_duplicates(subset=["lon", "lat"], keep="first").copy() + detectors["location"] = detectors.apply(lambda x: Point(x.lon, x.lat), axis=1) + detectors["geom"] = polygon + + intersect_dets = detectors[ + detectors.apply(lambda x: x.geom.contains(x.location), axis=1) + ]["detector_id"] + return forecast[forecast["detector_id"].isin(intersect_dets)] + + +def restrict_readings_to_network(network, readings, max_dist=5e-4): + + """Keep readings which lie close to network of interest""" + + detectors = readings.drop_duplicates(subset=["lon", "lat"], keep="first").copy() + + arr = detectors.apply( + lambda x: ox.distance.get_nearest_edge(network, (x.lat, x.lon), return_dist=True), + axis=1, + ) + + detectors["nearest_u"] = arr.apply(lambda x: x[0]) + detectors["nearest_v"] = arr.apply(lambda x: x[1]) + detectors["dist"] = arr.apply(lambda x: x[3]) + + detectors = detectors[detectors["dist"] <= max_dist] + detectors = detectors[ + ["detector_id", "lon", "lat", "nearest_u", "nearest_v", "dist"] + ] + + network_readings = readings.merge( + detectors, how="left", on=["detector_id", "lon", "lat"] + ) + return network_readings.dropna() + + +def aggregate_edge_event_count(network, network_df, non_zero_edges): + + """Aggregate sum of vehicle counts to each edge""" + + # Aggregates from detector level to grid=cell level + # stores in the graph and in the dataframe + # Update to be class G method + + agg_network = network.copy() + + edge_count_dict = {} + num_edges = 0 + + t_min = network_df["measurement_start_utc"].min() + t_max = network_df["measurement_end_utc"].max() + + t_ticks = pd.date_range(start=t_min, end=t_max, freq="H") + + nx.set_edge_attributes(agg_network, [], "counts") + nx.set_edge_attributes(agg_network, [], "baselines") + + for i, (source, target) in enumerate(non_zero_edges): + baselines = [] + counts = [] + for tick in t_ticks[:-1]: + edge_df = network_df[ + (network_df["nearest_u"] == source) + & (network_df["nearest_v"] == target) + & (network_df["measurement_start_utc"] == tick) + & (network_df["measurement_end_utc"] == tick + np.timedelta64(1, "h")) + ] + + edge_count = edge_df["count"].sum() / 1e6 + edge_baseline = edge_df["baseline"].sum() / 1e6 + + baselines.append(edge_baseline) + counts.append(edge_count) + + edge_count_dict[num_edges] = { + "source": source, + "target": target, + "measurement_start_utc": tick, + "measurement_end_utc": t_max, + "count": edge_count, + "baseline": edge_baseline, + } + + num_edges += 1 + + nx.set_edge_attributes(agg_network, {(source, target, 0): baselines}, name="baselines") + nx.set_edge_attributes(agg_network, {(source, target, 0): counts}, name="counts") + + if i % 50 == 0: + logging.info("Aggregation progress: %.2f%%", i * 100 / len(non_zero_edges)) + return agg_network, pd.DataFrame.from_dict(edge_count_dict, "index") + + +def path_event_count(network, path: list, t_min, t_max): + + """Sum vehicle counts on a given path""" + + baseline = 0 + count = 0 + + for i, _ in enumerate(path[:-1]): + + try: + edge = network[path[i]][path[i + 1]] + except KeyError: + edge = network[path[i + 1]][path[i]] + + baseline += sum(edge[0]["baselines"][t_min:t_max]) + count += sum(edge[0]["counts"][t_min:t_max]) + return baseline, count diff --git a/SpatialScan/network/path_generator.py b/SpatialScan/network/path_generator.py new file mode 100644 index 0000000..1e29d2e --- /dev/null +++ b/SpatialScan/network/path_generator.py @@ -0,0 +1,112 @@ +"""Generate paths on a osmnx network""" + +import numpy as np +import networkx as nx + + +class PathGenerator: + + """Base class for generating paths on a network using networkx/dijkstra""" + + def __init__( + self, + network: nx.MultiGraph, + min_path_length: float, + max_path_length: float, + detector_edges: list, + drop_duplicates=True, + drop_detectorless=True, + ): + self.network = network + self.min_path_length = min_path_length + self.max_path_length = max_path_length + self.detector_edges = detector_edges + self.drop_duplicates = drop_duplicates + self.drop_detectorless = drop_detectorless + self.network_paths = {} + + def generate_paths(self): + """Generate Paths""" + self.network_paths = self._find_all_paths() + if self.drop_duplicates: + self._drop_duplicate_paths() + if self.drop_detectorless: + self._drop_detectorless_paths() + + def _find_all_paths(self) -> dict: + """Find all possible paths with appropiate lengths on the network""" + + upper_paths = dict( + nx.all_pairs_dijkstra_path( + self.network, cutoff=self.max_path_length, weight="length" + ) + ) + lower_paths = dict( + nx.all_pairs_dijkstra_path( + self.network, cutoff=self.min_path_length, weight="length" + ) + ) + + useful_paths = upper_paths.copy() + + for source in upper_paths.copy().keys(): + for target in upper_paths[source].copy().keys(): + if ( + source in lower_paths.keys() + and target in lower_paths[source].keys() + ): + del useful_paths[source][target] + return useful_paths + + def _drop_detectorless_paths(self): + """Drop paths without detectors on them""" + + for source in self.network_paths.copy().keys(): + for target in self.network_paths[source].copy().keys(): + + path = self.network_paths[source][target] + path_graph = self.network.subgraph(path) + + if not any( + path_graph.has_edge(u, v) or path_graph.has_edge(v, u) + for (u, v) in self.detector_edges + ): + del self.network_paths[source][target] + + for source in self.network_paths.copy().keys(): + if not self.network_paths[source]: + del self.network_paths[source] + + def _drop_duplicate_paths(self): + """Drop duplicate paths ie. u -> v == v -> u""" + + seen_paths = set() + for source in self.network_paths.copy().keys(): + for target in self.network_paths[source].copy().keys(): + + if (source, target) in seen_paths: + del self.network_paths[target][source] + elif source == target: + del self.network_paths[target][source] + + seen_paths.add((source, target)) + seen_paths.add((target, source)) + + for source in self.network_paths.copy().keys(): + if not self.network_paths[source]: + del self.network_paths[source] + + def num_paths(self): + """Number of paths on the network""" + if len(self.network_paths) == 0: + raise TypeError("No Paths have been generated yet") + count = 0 + for source in self.network_paths.keys(): + count += len(self.network_paths[source]) + return count + + def random_path(self): + """Choose random path""" + source = int(np.random.choice(list(self.network_paths.keys()))) + target = int(np.random.choice(list(self.network_paths[source].keys()))) + return self.network_paths[source][target] From b5af9989e2879d6514945ce28580736710ee5875 Mon Sep 17 00:00:00 2001 From: Chance Haycock Date: Thu, 13 Aug 2020 11:30:21 +0100 Subject: [PATCH 3/5] black --- .../__pycache__/__init__.cpython-36.pyc | Bin 0 -> 151 bytes .../__pycache__/network_scan.cpython-36.pyc | Bin 0 -> 8212 bytes .../__pycache__/path_generator.cpython-36.pyc | Bin 0 -> 3235 bytes SpatialScan/network/network_scan.py | 43 +++++++++++++----- 4 files changed, 32 insertions(+), 11 deletions(-) create mode 100644 SpatialScan/network/__pycache__/__init__.cpython-36.pyc create mode 100644 SpatialScan/network/__pycache__/network_scan.cpython-36.pyc create mode 100644 SpatialScan/network/__pycache__/path_generator.cpython-36.pyc diff --git a/SpatialScan/network/__pycache__/__init__.cpython-36.pyc b/SpatialScan/network/__pycache__/__init__.cpython-36.pyc new file mode 100644 index 0000000000000000000000000000000000000000..94c326dc8cc75d05814c304bde4cb786f6688d77 GIT binary patch literal 151 zcmXr!<>do*&IfBe5OyS+yQfGEg z^>kPDWK~Nu?P@UFGA0N_a3DBsa3FvPUkD<&_{N39Ubqk(AcF%DMBvbo2!exy2>4!B zbx+UiYO#B!US?%xWo5p6?|VP>rAEWKxA9B9`#DYfcWvosp?)1#^4B0t^E9S2W1xGw zDh;5`G*}(3d270MTVrKb zd7`ljZ_P|^ot_^*c>#}q=4x<#6#@!&@zs(~q zf;1LrXbr*xJ_z^YnE7Z6;|@JD{}E^Cm;5RS*T4~txZ#;>^#@v4XC+pCQUzZs%*Ogk ztjZjeWmaQ#loi%sO_VlUVJ(!^A83~5Iau|YtQD)KMNeM3wY+LKjonTl#IYNX(|lKP z|4y3>AIGdGOdrczA>Zza#}&XJfl&z8fB_cRkYf;G&RaWtv1u9<+MK2`VG}*W>!}2 zT3P9jXjYovLVFX}%FNU?l&zVD5`4R$^;c(T@z-!y#@)KQ+hEjoPSH(8x70ipS^?9* zIeTyN((dk6{IehpuHL+V@6D^ep`+Sw1H4&M)s zJ8|dWD&vo?-Ur{qU_i=xs$Q_<>YOzR2JDWG#U|dmjw=DwYMNcv@n<-CORwvWZs<+J zM$4M9rrSnSZ&EM*CL4RX)X{5p=M3yXuhq{0B|Lf!SMp^L=o9poFt}$h^J!7LpiT8m z??b<5`czLXXt&-krg-g=HRnRQ?v+BCajC^#X2E;+H}>l;+I* zj>1Bq4XhftVoWaH9}nF}d_U|AIFxc6rHKnsR$SU^m*h62P4FJJrbON6k4VSWK(UG! zC@$oC`8nLS^?Y-u_a)a(W`S<06S#)a39VH#4E^*5KRS4^K@#x>E+r1z*lxPCJrati z=2Q>0DMVvx9{p=-0cb5&|GsfSo9VA>?|k?jV`@!HM}H1>lmPb{+Wur(rn`&)87t{sQ|=Bq8G}UXZ?sJy8Ns(S z0oT&%3VtX^ruLo@h4h!AnDL}tg&M?T(c!X;;Sfk^Jq!{qOG+>>Ne5^tEy8WjA?o); zFxnSY%pyp(C+A4r-I0hPT}kpa_ww!+yO%FNW7K!)GINtr(BWT5!yzvWPy%)A>~khd zVZ`{+v;W1fy(*0nlMaj`haVI5xY(vQTmq3*9*u|O@o=wEsJa5+L|O;@IFT*2uCwwi zO#rqax%6ykdcJIuIMGy-z7G3)O9VqMtIC`(m{CER)~+kbmzC5{!_GnC)r7(upGFi+ z$hl+Dt|VyZKoFwg5yitqe1+cg%S65k;yEW?+8Lk#>Vti+W01l`&D=07-7>a;wM`&~ zp_Vxjl?#);4y6Wha(1>5bgFrjt4vh97?mPz_>c#Wwcw&rB=-rwzIlsKX z=OwMz9(NDdvYSZe` z`)VC)#X7{4P=w~J0%e>7{kjG`!F^VPajR$b{)+gw!yjdJ=HQ)hDW?s37qA8>k~LWE zgAx#B2}kl>(HLfSF(QP^iS=P0tp1dC_@f6}R%VS&jhP3AcsE_`udyatR)96@m=z-n zw0Re3qxCoX=du!ODTLWnPBmOzh{F2h({q23zyxZ&kb{LO5ijA8+NStL6oNc=@g*V( zTLf9WiK7I^Prv@%sS=z?`Mc9AE^ZWNd<{d|C8aL8mMoA&s4iVv5tLF;2r#*Hqp(?s zh&}=_;uU)F*NKn`Qlz~lzCrctL>365?3=VY@#t6@(UDw-+0Tt&kODCuj+MEG?>Y=p zaf2RJCY{7x`~rwH;blv^6Gv%)m_?TOaFiZPi&!KUkRnKn#Z4l&hC(QBcKu_}AV5cJ{-lSQujfi>Vyug`MJUa_1OP0Ry zF3As8wmQWO3!rff%UsjX18}yD7xW#296ATCiec0XumSiOHu_z&uD@Wk^vR3MdVFHt z&%OroQO29pAqvv!c8y$WI#{^Wj_(g+Hbyw$%ewC$j)OtbQ}KNkcYI$tpvs35S81NN zh}+ARfWe=k4kAFMpQpO=z$ktsz)RLZun7yHo8G)+J9gRryj>F1uU7e! zpE)4wC{)AqW^pOdqT?sf>Ko)%Bpa#8P`(cz7ET9hQzJ9_Wu|;xsspGj0EBh``07_E zrld#`YFL?k`W8<+`)&~fxkaSoQn*Ai3T>FI-$Fd)vY01s6sImff`1zR2R7W+J&TVL zL|1z^?%Y(0>DfIVlXGy~wxtm#N-bqMOeiLk&Tj;RF~2E9ETq|k6|PeFq~bNN48iG$ zM0^u_!=aMX<~53CCHau*AP-+5jEjTFPhbV3{}Na31DETF793bJLG26NCl->dC>Bz< z!Ql1-PYrlbdTs-li9j5F+EZ(4qy{iD4+P4S-JDH@N+}B#0w?ki65Ja{EOE)WyKyp% zjwYY~tsSyGkAvf#-`*LuOYis^h37mFn(r+ z`KX3NM0}^Dy)pUh&Rd+1PVyIz5yR}egOHy&QFwGQLTY3)WT_j9B!x>t5<_kvA1x$` zOc6=deOw8-6(^z7j8bm;9weDGSoDb;fWRNtkIsm%v|kwn!-p)m_V*Z(;L;|$6h~6N zEc3q+QJ(P#kk6((VX4wZ2TsiQi7BJr9C zdXY;&YJ@T*$(c7jC)H9g;_qk)#CXlx)Hk8^@L%DvipOaA@&EkrKmYO0!|UJZ&=Jv( zl%I?}-6wwmf>{5f5(jzd8;hK!olD)+$}GrXmE;VLnyF-uT9o9emO7-7S{jQnO|)P( z)&%w%2arKXrHKqDAP1S;I&1cT!f4Y;Mp-NG=^?0=pwMx_QYTRigU|^Zv_XUdrkGRP#O~S^fX! zmGgO}MMoisl}ID;E?`TV`|)F-dxX!`o8L75m_&a~D;w1lg1^jC;{Y_?g%SARNf!qDgxJt4Y+v&9fs zWwit2I;G~Zo~&8av7W4%+WjgZqCPXW;A&mhpxi6j3Tt4k&P-#?i`uk>`<12pRm5nu zY?YGsSp&6}s!=Yei72d!ab@J-2Y6mU3*noE^&dLzDVRG$|GM${2CF8+z{W5 zyc*?Y7$2m;mtRKgHXb6EI~pMo_S|5ty_V~lGxy9`*K0=jC^O-r2%W?OOz*7_EAy;r zo)4;09H(NdrFio;E$7=rPWbUY-Tp3-H$l9Lnh39PPGKxRlc*8Oj4{itSWHp}ak1PL#BXE7U=oU1E`9 t7V!zbcAgQq23~2&iy#!9*!a$r|55A+NjF=+qit!rQ@87O)7iGa`hQNV^b7z1 literal 0 HcmV?d00001 diff --git a/SpatialScan/network/__pycache__/path_generator.cpython-36.pyc b/SpatialScan/network/__pycache__/path_generator.cpython-36.pyc new file mode 100644 index 0000000000000000000000000000000000000000..a4a23b1336e3374af5be29510ac0102e51d095bb GIT binary patch literal 3235 zcma)8&u`<#6`mQAGA;jH$MM=}8Zhl*5v7QoA}xvnYa>Y9-J(d10PA)OVH1|C8CjHR zQl24qWkI^E475G;QskIJuf6uq;kBm%{Rf&;zc-|0%iRK%z&A6G!};;P?|n0$tgSVF z`R|{E|5|74AMDo8LH$!a<{^w^Z6>)C32zHqI&Jp}la6$+nRKV3?Md;JwY<^$PgJUO zl&fG6Pa@UF@=xsGy8ZpclP>++u5Lb)N| z!V=#u*S^Hc>Sbnnbo^OlRM1HxV}fp`gR^P_ae8)p3-#iI!qAJlH9U}U|J>v{I&Hb- zQtd4_F5)z#?T3j<&+;A|(Qv`BE>*5Nd8R`p&y?vbVc+zlQsDUA$2wgWid$N4!|+xy~eb zYxVl=N=<92bd5^7ZyHY1I<*h;H>OumeD|w^-x#IML8lj`9o37jI$7uZK&qDqPZ5YX zN}hJ2bfI2Pc3|UF$2`;fgR9aH!#IueFr<+UXvx@egWvo|I|~5r&$@Rwb;Xa$c@#fH zW6b(&N>M^I`$F%+I#YB*>FS5D9k&)DVTXKd^sq6`%F;mS4Iu z?&PINhvoeWzgU}ambNy~6~k?=P}Y=gCj$%Jb9H4}4GUCzr5Fs9u3%}eB-yKZW3@1< zThzQL8zw7sN0p5{(q}5a1F$A8Kfq(&fnko|F5jxG>x(VEBUbq)%F*{1QNF{zIk>c1 z0+$lr{1s+60I00b0f?|LeLdmB-;r<30mF^|vk6haddFf9@S^4P-7%=}CzgAEU|miW zxc1?19J48(a`VB>7Xa`B0J!;#4XyXq$U76tF*~w(UDeS0^)s9mfZ82BI0lxruoAji zkq0v|oZ)&g-y!XccKHb&Q)Tj3KX@)Cd?Nm6XBmsB7zDNss!qjkk;j)}riHxpdy%mk z+j2^=|8dESlD|ZrGzTQlXF3}6ZW8djBDDz8~mS`8(;`)7ZlSt@I zm_s~f2ZlM@d{?aVRUvp2c3W(t96kL1Ry+?gg_!tnRpRfc-~F1OyUerwa&>_g~l^4l-RFvUHo?oD2Q?+g@l1bb`K^a>y{XFhbS3+y}LG zYHHn{eaDj4&TYP_Rsm{#18P;m-Lll0pK&$GVznO>!AD2I%i!=32DRd(<&};>elCIo zG^DJ7Xg7^g)hTZD9F)5Wy@6PLiO1~1FsH!@j#a)6v0eo%I~5#1d?OsQv$$f<91EiS z9E)Q#ayq7vfRZDMJGqdE*=w|b&jWf>V?)t_idS*jL(2XP2wpyf#!_^g(R(Mw#S5i_ ztb22I*;*nb#C>Y9M_`0 zM`sa<)y;`Edt(Lg<}JyeM?PXQ;=@Nn4?9sV94?Je&5R?m}`I6ZA`(v>KUE>swn z%`m*kWRXyPH4HC{D5>8pg`vzkVW{t8;IjGaBFWpF literal 0 HcmV?d00001 diff --git a/SpatialScan/network/network_scan.py b/SpatialScan/network/network_scan.py index 8552fe9..1080ec9 100644 --- a/SpatialScan/network/network_scan.py +++ b/SpatialScan/network/network_scan.py @@ -16,6 +16,7 @@ class NetworkScan(PathGenerator): """Functionality to carry out Scan on Road Networks""" + def __init__( self, forecast, @@ -39,7 +40,11 @@ def __init__( super().__init__( self.network, min_path_length, max_path_length, detector_edges, **kwargs ) - logging.info("Calculating paths of lengths between %f and %f", self.min_path_length, self.max_path_length) + logging.info( + "Calculating paths of lengths between %f and %f", + self.min_path_length, + self.max_path_length, + ) self.generate_paths() self.t_min = forecast["measurement_start_utc"].min() @@ -113,7 +118,6 @@ def highest_path(self): """Path with highest EBP score""" return self.results.iloc[0]["path"] - def aggregate_results_to_edges(self): """Calculate mean EBP score per edge""" @@ -130,8 +134,15 @@ def aggregate_results_to_edges(self): agg_dict = {} for i, (source, target) in enumerate(self.detector_edges): - spatial_cond = (self.results["path"].astype(str).str.contains("{}, {}".format(source, target)))\ - | (self.results["path"].astype(str).str.contains("{}, {}".format(target, source))) + spatial_cond = ( + self.results["path"] + .astype(str) + .str.contains("{}, {}".format(source, target)) + ) | ( + self.results["path"] + .astype(str) + .str.contains("{}, {}".format(target, source)) + ) for tick in t_ticks[:-1]: @@ -153,11 +164,15 @@ def aggregate_results_to_edges(self): num_aggs += 1 nx.set_edge_attributes( - agg_network, {(source, target, 0): score}, name="score_t={}".format(tick) + agg_network, + {(source, target, 0): score}, + name="score_t={}".format(tick), ) if i % 50 == 0: - logging.info("Aggregation progress: %.2f%%", i * 100 / len(self.detector_edges)) + logging.info( + "Aggregation progress: %.2f%%", i * 100 / len(self.detector_edges) + ) self.agg_results = pd.DataFrame.from_dict(agg_dict, "index") self.network = agg_network @@ -170,7 +185,7 @@ def get_borough_polygon( """Fetch borough polygon information""" if not os.path.isfile(boroughs_file): - raise ValueError('File does not exist') + raise ValueError("File does not exist") london_df = gpd.read_file(boroughs_file) london_df = london_df.to_crs(epsg=4326) @@ -186,7 +201,7 @@ def get_network_from_polygon(polygon): network = ox.graph_from_polygon( polygon, network_type="drive", simplify=True, custom_filter=roi ) - return nx.MultiGraph(network) # Note - not directional + return nx.MultiGraph(network) # Note - not directional def restrict_readings_to_polygon(forecast, polygon): @@ -210,7 +225,9 @@ def restrict_readings_to_network(network, readings, max_dist=5e-4): detectors = readings.drop_duplicates(subset=["lon", "lat"], keep="first").copy() arr = detectors.apply( - lambda x: ox.distance.get_nearest_edge(network, (x.lat, x.lon), return_dist=True), + lambda x: ox.distance.get_nearest_edge( + network, (x.lat, x.lon), return_dist=True + ), axis=1, ) @@ -278,8 +295,12 @@ def aggregate_edge_event_count(network, network_df, non_zero_edges): num_edges += 1 - nx.set_edge_attributes(agg_network, {(source, target, 0): baselines}, name="baselines") - nx.set_edge_attributes(agg_network, {(source, target, 0): counts}, name="counts") + nx.set_edge_attributes( + agg_network, {(source, target, 0): baselines}, name="baselines" + ) + nx.set_edge_attributes( + agg_network, {(source, target, 0): counts}, name="counts" + ) if i % 50 == 0: logging.info("Aggregation progress: %.2f%%", i * 100 / len(non_zero_edges)) From a75be562a0796afccb43ae2a1b3d8225ba4d46bd Mon Sep 17 00:00:00 2001 From: Chance Haycock Date: Thu, 13 Aug 2020 11:31:06 +0100 Subject: [PATCH 4/5] remove useless pycache --- .../network/__pycache__/__init__.cpython-36.pyc | Bin 151 -> 0 bytes .../__pycache__/network_scan.cpython-36.pyc | Bin 8212 -> 0 bytes .../__pycache__/path_generator.cpython-36.pyc | Bin 3235 -> 0 bytes 3 files changed, 0 insertions(+), 0 deletions(-) delete mode 100644 SpatialScan/network/__pycache__/__init__.cpython-36.pyc delete mode 100644 SpatialScan/network/__pycache__/network_scan.cpython-36.pyc delete mode 100644 SpatialScan/network/__pycache__/path_generator.cpython-36.pyc diff --git a/SpatialScan/network/__pycache__/__init__.cpython-36.pyc b/SpatialScan/network/__pycache__/__init__.cpython-36.pyc deleted file mode 100644 index 94c326dc8cc75d05814c304bde4cb786f6688d77..0000000000000000000000000000000000000000 GIT binary patch literal 0 HcmV?d00001 literal 151 zcmXr!<>do*&IfBe5OyS+yQfGEg z^>kPDWK~Nu?P@UFGA0N_a3DBsa3FvPUkD<&_{N39Ubqk(AcF%DMBvbo2!exy2>4!B zbx+UiYO#B!US?%xWo5p6?|VP>rAEWKxA9B9`#DYfcWvosp?)1#^4B0t^E9S2W1xGw zDh;5`G*}(3d270MTVrKb zd7`ljZ_P|^ot_^*c>#}q=4x<#6#@!&@zs(~q zf;1LrXbr*xJ_z^YnE7Z6;|@JD{}E^Cm;5RS*T4~txZ#;>^#@v4XC+pCQUzZs%*Ogk ztjZjeWmaQ#loi%sO_VlUVJ(!^A83~5Iau|YtQD)KMNeM3wY+LKjonTl#IYNX(|lKP z|4y3>AIGdGOdrczA>Zza#}&XJfl&z8fB_cRkYf;G&RaWtv1u9<+MK2`VG}*W>!}2 zT3P9jXjYovLVFX}%FNU?l&zVD5`4R$^;c(T@z-!y#@)KQ+hEjoPSH(8x70ipS^?9* zIeTyN((dk6{IehpuHL+V@6D^ep`+Sw1H4&M)s zJ8|dWD&vo?-Ur{qU_i=xs$Q_<>YOzR2JDWG#U|dmjw=DwYMNcv@n<-CORwvWZs<+J zM$4M9rrSnSZ&EM*CL4RX)X{5p=M3yXuhq{0B|Lf!SMp^L=o9poFt}$h^J!7LpiT8m z??b<5`czLXXt&-krg-g=HRnRQ?v+BCajC^#X2E;+H}>l;+I* zj>1Bq4XhftVoWaH9}nF}d_U|AIFxc6rHKnsR$SU^m*h62P4FJJrbON6k4VSWK(UG! zC@$oC`8nLS^?Y-u_a)a(W`S<06S#)a39VH#4E^*5KRS4^K@#x>E+r1z*lxPCJrati z=2Q>0DMVvx9{p=-0cb5&|GsfSo9VA>?|k?jV`@!HM}H1>lmPb{+Wur(rn`&)87t{sQ|=Bq8G}UXZ?sJy8Ns(S z0oT&%3VtX^ruLo@h4h!AnDL}tg&M?T(c!X;;Sfk^Jq!{qOG+>>Ne5^tEy8WjA?o); zFxnSY%pyp(C+A4r-I0hPT}kpa_ww!+yO%FNW7K!)GINtr(BWT5!yzvWPy%)A>~khd zVZ`{+v;W1fy(*0nlMaj`haVI5xY(vQTmq3*9*u|O@o=wEsJa5+L|O;@IFT*2uCwwi zO#rqax%6ykdcJIuIMGy-z7G3)O9VqMtIC`(m{CER)~+kbmzC5{!_GnC)r7(upGFi+ z$hl+Dt|VyZKoFwg5yitqe1+cg%S65k;yEW?+8Lk#>Vti+W01l`&D=07-7>a;wM`&~ zp_Vxjl?#);4y6Wha(1>5bgFrjt4vh97?mPz_>c#Wwcw&rB=-rwzIlsKX z=OwMz9(NDdvYSZe` z`)VC)#X7{4P=w~J0%e>7{kjG`!F^VPajR$b{)+gw!yjdJ=HQ)hDW?s37qA8>k~LWE zgAx#B2}kl>(HLfSF(QP^iS=P0tp1dC_@f6}R%VS&jhP3AcsE_`udyatR)96@m=z-n zw0Re3qxCoX=du!ODTLWnPBmOzh{F2h({q23zyxZ&kb{LO5ijA8+NStL6oNc=@g*V( zTLf9WiK7I^Prv@%sS=z?`Mc9AE^ZWNd<{d|C8aL8mMoA&s4iVv5tLF;2r#*Hqp(?s zh&}=_;uU)F*NKn`Qlz~lzCrctL>365?3=VY@#t6@(UDw-+0Tt&kODCuj+MEG?>Y=p zaf2RJCY{7x`~rwH;blv^6Gv%)m_?TOaFiZPi&!KUkRnKn#Z4l&hC(QBcKu_}AV5cJ{-lSQujfi>Vyug`MJUa_1OP0Ry zF3As8wmQWO3!rff%UsjX18}yD7xW#296ATCiec0XumSiOHu_z&uD@Wk^vR3MdVFHt z&%OroQO29pAqvv!c8y$WI#{^Wj_(g+Hbyw$%ewC$j)OtbQ}KNkcYI$tpvs35S81NN zh}+ARfWe=k4kAFMpQpO=z$ktsz)RLZun7yHo8G)+J9gRryj>F1uU7e! zpE)4wC{)AqW^pOdqT?sf>Ko)%Bpa#8P`(cz7ET9hQzJ9_Wu|;xsspGj0EBh``07_E zrld#`YFL?k`W8<+`)&~fxkaSoQn*Ai3T>FI-$Fd)vY01s6sImff`1zR2R7W+J&TVL zL|1z^?%Y(0>DfIVlXGy~wxtm#N-bqMOeiLk&Tj;RF~2E9ETq|k6|PeFq~bNN48iG$ zM0^u_!=aMX<~53CCHau*AP-+5jEjTFPhbV3{}Na31DETF793bJLG26NCl->dC>Bz< z!Ql1-PYrlbdTs-li9j5F+EZ(4qy{iD4+P4S-JDH@N+}B#0w?ki65Ja{EOE)WyKyp% zjwYY~tsSyGkAvf#-`*LuOYis^h37mFn(r+ z`KX3NM0}^Dy)pUh&Rd+1PVyIz5yR}egOHy&QFwGQLTY3)WT_j9B!x>t5<_kvA1x$` zOc6=deOw8-6(^z7j8bm;9weDGSoDb;fWRNtkIsm%v|kwn!-p)m_V*Z(;L;|$6h~6N zEc3q+QJ(P#kk6((VX4wZ2TsiQi7BJr9C zdXY;&YJ@T*$(c7jC)H9g;_qk)#CXlx)Hk8^@L%DvipOaA@&EkrKmYO0!|UJZ&=Jv( zl%I?}-6wwmf>{5f5(jzd8;hK!olD)+$}GrXmE;VLnyF-uT9o9emO7-7S{jQnO|)P( z)&%w%2arKXrHKqDAP1S;I&1cT!f4Y;Mp-NG=^?0=pwMx_QYTRigU|^Zv_XUdrkGRP#O~S^fX! zmGgO}MMoisl}ID;E?`TV`|)F-dxX!`o8L75m_&a~D;w1lg1^jC;{Y_?g%SARNf!qDgxJt4Y+v&9fs zWwit2I;G~Zo~&8av7W4%+WjgZqCPXW;A&mhpxi6j3Tt4k&P-#?i`uk>`<12pRm5nu zY?YGsSp&6}s!=Yei72d!ab@J-2Y6mU3*noE^&dLzDVRG$|GM${2CF8+z{W5 zyc*?Y7$2m;mtRKgHXb6EI~pMo_S|5ty_V~lGxy9`*K0=jC^O-r2%W?OOz*7_EAy;r zo)4;09H(NdrFio;E$7=rPWbUY-Tp3-H$l9Lnh39PPGKxRlc*8Oj4{itSWHp}ak1PL#BXE7U=oU1E`9 t7V!zbcAgQq23~2&iy#!9*!a$r|55A+NjF=+qit!rQ@87O)7iGa`hQNV^b7z1 diff --git a/SpatialScan/network/__pycache__/path_generator.cpython-36.pyc b/SpatialScan/network/__pycache__/path_generator.cpython-36.pyc deleted file mode 100644 index a4a23b1336e3374af5be29510ac0102e51d095bb..0000000000000000000000000000000000000000 GIT binary patch literal 0 HcmV?d00001 literal 3235 zcma)8&u`<#6`mQAGA;jH$MM=}8Zhl*5v7QoA}xvnYa>Y9-J(d10PA)OVH1|C8CjHR zQl24qWkI^E475G;QskIJuf6uq;kBm%{Rf&;zc-|0%iRK%z&A6G!};;P?|n0$tgSVF z`R|{E|5|74AMDo8LH$!a<{^w^Z6>)C32zHqI&Jp}la6$+nRKV3?Md;JwY<^$PgJUO zl&fG6Pa@UF@=xsGy8ZpclP>++u5Lb)N| z!V=#u*S^Hc>Sbnnbo^OlRM1HxV}fp`gR^P_ae8)p3-#iI!qAJlH9U}U|J>v{I&Hb- zQtd4_F5)z#?T3j<&+;A|(Qv`BE>*5Nd8R`p&y?vbVc+zlQsDUA$2wgWid$N4!|+xy~eb zYxVl=N=<92bd5^7ZyHY1I<*h;H>OumeD|w^-x#IML8lj`9o37jI$7uZK&qDqPZ5YX zN}hJ2bfI2Pc3|UF$2`;fgR9aH!#IueFr<+UXvx@egWvo|I|~5r&$@Rwb;Xa$c@#fH zW6b(&N>M^I`$F%+I#YB*>FS5D9k&)DVTXKd^sq6`%F;mS4Iu z?&PINhvoeWzgU}ambNy~6~k?=P}Y=gCj$%Jb9H4}4GUCzr5Fs9u3%}eB-yKZW3@1< zThzQL8zw7sN0p5{(q}5a1F$A8Kfq(&fnko|F5jxG>x(VEBUbq)%F*{1QNF{zIk>c1 z0+$lr{1s+60I00b0f?|LeLdmB-;r<30mF^|vk6haddFf9@S^4P-7%=}CzgAEU|miW zxc1?19J48(a`VB>7Xa`B0J!;#4XyXq$U76tF*~w(UDeS0^)s9mfZ82BI0lxruoAji zkq0v|oZ)&g-y!XccKHb&Q)Tj3KX@)Cd?Nm6XBmsB7zDNss!qjkk;j)}riHxpdy%mk z+j2^=|8dESlD|ZrGzTQlXF3}6ZW8djBDDz8~mS`8(;`)7ZlSt@I zm_s~f2ZlM@d{?aVRUvp2c3W(t96kL1Ry+?gg_!tnRpRfc-~F1OyUerwa&>_g~l^4l-RFvUHo?oD2Q?+g@l1bb`K^a>y{XFhbS3+y}LG zYHHn{eaDj4&TYP_Rsm{#18P;m-Lll0pK&$GVznO>!AD2I%i!=32DRd(<&};>elCIo zG^DJ7Xg7^g)hTZD9F)5Wy@6PLiO1~1FsH!@j#a)6v0eo%I~5#1d?OsQv$$f<91EiS z9E)Q#ayq7vfRZDMJGqdE*=w|b&jWf>V?)t_idS*jL(2XP2wpyf#!_^g(R(Mw#S5i_ ztb22I*;*nb#C>Y9M_`0 zM`sa<)y;`Edt(Lg<}JyeM?PXQ;=@Nn4?9sV94?Je&5R?m}`I6ZA`(v>KUE>swn z%`m*kWRXyPH4HC{D5>8pg`vzkVW{t8;IjGaBFWpF From 9e901f37e570a6d3115a59d67cb58721fe29f21e Mon Sep 17 00:00:00 2001 From: Chance Haycock Date: Thu, 13 Aug 2020 11:32:09 +0100 Subject: [PATCH 5/5] Update gitignore --- .gitignore | 1 + 1 file changed, 1 insertion(+) diff --git a/.gitignore b/.gitignore index b676934..5f772aa 100644 --- a/.gitignore +++ b/.gitignore @@ -2,6 +2,7 @@ .vscode* *__pycache__/* */__pycache__/* +network/__pycache__/* .ipynb* region_event_try.py *.csv