{'fraction_non_alphanumeric': 0.028607144574258378, 'fraction_numerical': 0.010865907014760324, 'mean_word_length': 4.639013148976549, 'pattern_counts': {'":': 0, '<': 2, '<?xml version=': 0, '>': 3, 'https://': 1, 'lorem ipsum': 0, 'www.': 1, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Intertemporal Community Detection in Bikeshare Networks

We investigate the changes in the patterns of usage in the Divvy bikeshare system in Chicago from 2016-2018. We devise a community detection method that finds clusters of nodes that are increasing, decreasing, or stable in connectivity across time. We use an iterative testing approach that is augmented by trend testing and a novel temporal false-discovery-rate correction. Results show stark geographical patterns in clusters that are growing and declining in relative bike-share usage across time and may elucidate latent economic or demographic trends.


Introduction
Many of the community detection techniques in networks have been applied to static networks [Girvan and Newman, 2002]. However, many of networks exhibit dynamic properties. In cities, these networks include commuting patterns over time [Patuelli et al., 2010], location based social networks [Assem et al., 2016], taxicab travel patterns [Liu et al., 2015] and cell phone call records [Reades et al., 2009]. Understanding the structure of these networks and how it varies over time is important to understand the rhythms of the city and manage its infrastructure.
Bikeshare systems have become popular in many cities (about 2000 cities in 2019) across the world as a means of providing mobility solutions. In 2018, according to National Association of City Transportation Officials, 36.5 million trips were completed in over 100 cities in the United States. Many of these systems have stations, where users can rent the bikes and deposit them at another station at the end of the trip. In recent years, dockless systems (without fixed stations) are becoming more prevalent. Nonetheless, the nature of these systems allow for the system operator to track the precise origins and destinations of the trip by the time of day and day of the week. These datasets are rich and are amenable to network analysis because of their structure and is remarkably comprehensive to track mobility patterns within a city, despite the limitations of use based on age, social status and uneven distribution of stations as well maintenance schedules.
In this paper, we use datasets from Divvy Bikeshare in Chicago that contains fine-grained temporal data for all linkages between all bike stations for several years from its inception in 2013. Analysis of time-varying weighted graphs may allow us to gain more insight to the nature of the city as a complex process of microcosmic activity flows. Existing time-series as well as network -based analysis of bikeshare data exists [Austwick et al., 2013, Cazabet et al., 2017, but they do not fully take into account the dependencies induced by the network structures and temporal trends. In this work, we model the system as a time-varying weighted network. Bikeshare stations are represented as nodes and the total number of trips between the two nodes are represented as weights on the edges. We Figure 1: Communities of stations detected using a weighted configuration model extraction method (CCME) for aggreated number of trips between two stations for the whole dataset.
seek to find clusters of significantly connected stations across time with connectivities that may be increasing, decreasing, or stable across time. Identification of such clusters allow us to understand the nature of geographical, economic and cultural relationships within the city. However, clusters that are increasing or decreasing in connectivity across time may reveal latent signals structural changes within these clusters.
A naive way to approach the problem of community detection on a time-series of graphs is to aggregate the edge weight in the graphs and apply community detection on the graph, ignoring the temporal variation. Such an approach may find communities, but does not capture the trends or changes. Groups of nodes (stations) with few riders amongst themselves at one time but have grown more connected will be on par with clusters of nodes that have stayed relatively constant. For example, applying community detection procedure developed by Palowitch et al. [2018] on the aggregation, we find three broad clusters The method is able to split up the stations across very general geographical partitions ( Figure 1) with minimal overlaps. However, even after adjusting the tuning parameters, the method is unable to detect fine grained variation that map to neighborhoods (except for the downtown Loop).

Hypothesis and Data
We use the open source dataset to explore trips between 2016-2018 on the bikeshare system. Overall ridership declined from 3.6 million trips in 2017 to 3.39 million in 2018 ( Figure 2). The number of riders in 2018 is slightly lower number than ridership in 2017 at 3.4 million. Such a decline may be indicative of overarching trends in demographic shifts that systematically shift ridership, or auxiliary factors affecting ridership such as weather and construction projects, or a combination of the two. Overarching trends may be attributed to micro-level migration patterns in particular neighborhoods, for example a decrease in ridership may be due to instability of residency among residents. We aim to find clusters of stations that show broader demographic patterns that are consistent across clusters.
A further question that arises is just how global or local are trends in ridership. In particular, is the 2018 decline equally distributed across all stations, or are there certain areas that declined more than others, while certain other clusters have actually increased? We aim to answer both of these question by developing a method for intertemporal community detection that finds groups of nodes that: 1. are consistently connected across time.
2. exhibit a positive, negative, or neutral trend in connectivity; in particular, we hope to find decreasing clusters that may elucidate the overall decline in 2018.
To do so, we use the assumptions of the weighted configuration model as posited in Palowitch et al. [2018] and He et al. [2019]. We use normal approximation assumption made in the weighted configuration model to scale and to a certain extent whiten the time-series of node-set connectivities. By centering and scaling each instance of relative connectivity, we remove much of the overall graph effects that represent network-wide signals at a given time-point, such as weather and other city-wide phenomena at a given time-point. Though normalizing the node-set connectivties remove some of the seasonality, some autoregressive behavior still exist: we ignore these effects and simply fit a trend, but such behavior should be accounted for in future work.
Data across 547 Divvy stations that are active from June 2016 -June 2018 are used for this analysis. We use data from June 2016 because stations in the neighboring cities Evanston and Oak Park were only opened in June 2016, and Oak Park stations were closed in mid-2018. As such, the chosen time period accounts for perhaps the greatest number of stations across all time points. The publicly released datasets include trip start and stop times as well as stations and minimal attributes of the rider. We aggregate the trips between two stations by week and use them as edge weights between the two nodes.

Methods
In this section we describe a method based on iterative testing of node-set connectivities to extract statistically significant communities cross time [He et al., 2019, Palowitch et al., 2018. Such an approach is coupled with time-decay adjustments and with equivalence testing for trends [Schuirmann, 1987, Dixon andPechmann, 2008]: positive, negative or neutral trends are pre-specified categories before clustering with a bounding parameter U that represents the energy barrier of the slope of the trend. We introduce an iterative testing framework that finds clusters that are both significantly connected across time but also exhibiting increasing, decreasing, and stable trends in connectivity.

Configuration Null Model
In Palowitch et al. [2018], the weights on the edges incident on each node u are modeled, under the null hypothesis of connectivity, as where each du represents the total strengths from node u, where each su represents the total strengths from node u, and sT , dT represents the global sum of strengths and degrees.
Communities are sets of vertices that have edges that are significantly interconnected within the set but not connected outside the set. Palowitch et al. [2018], He et al. [2019] use the above null model to identify communities within a single graph. We now define a null model for the trend of normalized weights Zt(B), which are the sums of weights in a given node set B on random graphs Gt with fixed degrees du and strengths su and observed graphs Gt whose expected degrees and strengths are du and su [He et al., 2019]. We start with the definition of a node-set connectivity measure as defined in Palowitch et al. [2018] with means and variances detailed derivations of these values can be found in the text of Palowitch et al. [2018].

Normalizing Edges
One consequence of modeling network time-series with the weighted configuration model paradigm is that one is able to place edge-weights in a relative scale when they are normalized by their expectations and variances. That we can scale edge weights is especially useful in the time-series setting. Seasonal effects dominate much of the variation in the Divvy data. However, normalizing edge weights removes much of the global graph effects at given time t. Empirically, such a procedure removes much of the seasonality. The κ parameter is shown to be highly seasonal ( Figure 3). Normalization in effect takes the overall graph effects and residualizes for the effects of community interactions across time, as seen by the high seasonality in the variance terms. Effects of this normalization may be seen in a later section in Figure 6.

Intertemporal Null Model
The random variables ξuv are constructed as to satisfy the conditions of the weighted configuration model for the goal of community detection. However, because the domain of analysis for these previous studies is a single weighted graph, there is limited amount of information available to draw detailed conclusions or interpretations about the parameters. When presented with a set or time series of observations for each Wuv, we can learn more about the time-varying behaviors of network structure.
We significantly extend the model used in Palowitch et al. [2018], He et al. [2019] by developing an intertemporal null model wherein a baseline model is extracted from a time-series of registered networks and such a model iteratively subject to hypothesis tests for local deviance. After defining these key components of an intertemporal null model, we detect communities (clusters of stations) across the time-series of bikeshare networks that are significant compared to a null model if the trend component is significantly different than that of the baseline model. These communities signal regions of significant signals across time.
We now define the intertemporal null model for a given B that is significantly defined as a community (as in Palowitch et al. [2018], He et al. [2019]) in all time-points contingent on whether we assume the nodes are increasing, decreasing, or stable in connectivity. The notion of significance at all time points is explained in the following section. We search for communities that are not necessarily significantly connected at all time points, but which are • if increasing, significantly connected at time t = 1, but not necessarily as t grows larger.
• if decreasing, significantly connected at time t = T , but not necessarily when t is early, • if stable (or neutral), significantly connected at all time points.
Within a group of nodes B, we posit that a time-series of relative connectivity may be decomposed into trend and variation components. Trend denotes the presence of a significantly increasing or decreasing time-trend in the relative connectivity between all nodes in set B. Variation in this context denotes components of the node-set connectivity that do not vary systematically across time, i.e. with a constant mean.

Null Model for Node-Set Connectivity
We define an intertemporal null model for significance of connectivity in (4). This model is a simple extension of (1), which was the null model for a single graph.
To search for intertemporal clusters, we investigate not only the significance of connectivities of Bt to vt , but also their connectivities across time t. We define the intertemporal normalized connectivity score Zt(v, B) as: Finally, the sum of node-scores represent the overall significance of a cluster. Note that such an agglomeration no longer provides an adequate normal approximation, but does show the total connectivity amongst nodes in a cluster [Palowitch et al., 2018]. We define the sum of node-scores as Z(B) to gauge the significance of the time-trend of a cluster:

Null Model for Trend
For a given community B that is significantly connected across time t = 1, ..., T , we write the vector of node-set connectivity as Z(B), which spans from t = 1, ..., T , as where βu,Bt represents the trend component which is linearly dependent on time and V(B) represents the variation component that is stationary across time, but may have autoregressive or movingaverage characteristics.
The null condition for a set B is that it yields constant connectivity from all of its members through time, but in the case of a community that increases or decreases in relative connectivity. The null hypothesis for all v ∈ B is as follows: To detect clusters of nodes that exhibit increasing or decreasing connectivity, we utilize an iterative testing approach that, contingent on the hypothesis, finds clusters that are significantly connected and are increasing, decreasing, or neutral in connectivity across time. We split the community detection procedure into three parts corresponding to pre-specified trend hypotheses and perform the algorithm three times with procedures that are tailored to each hypothesized trend behavior of clusters.
Conversely, if the hypothesis in question is whether the trend is negligible, then the hypothesis in regard to the trend is a two sided test.
H0 :βv,B = 0 H1 :βv,B = 0 However, because the hypothesis in question is whether the time trend is zero, such a test must be inverted in order to capture the significance of the anomalous setting in which the trend is equal to zero. Such a test may be very hard to approach directly, so we therefore invoke the methods used in Dixon and Pechmann [2008] to approximate such a test using a pre-selected interval about zero. In this application, we use a symmetric interval [−U, U ].

Identifying Nodes that are Significantly Bordering Across Time
We use iterative testing to identify nodes that are significantly connected to its neighbors through time. Methods developed in previous literature [Palowitch et al., 2018, He et al., 2019 have applied this method to a fixed graph G. We use the same method of deriving a node-set p-value as in Palowitch et al. [2018], which describes the likelihood that v ∈ B is significantly connected to u for a fixed node u and a given set B.
We extend this methodology for a time-series of graphs Gt. At a given iteration k > 1, For fixed time t, with a node-set B k,t and a bordering node u, we determine the p-value of node-set connectivity by computing the connectivity measure from (2): and determine a p-value for each v ∈ B k,t by normalizing the test statistic as in (6) p(u, B k,t , Gt) = P(S(u, B k,t , Gt) > S(u, B k,t , Gt)).
For each time point, the p-value associated with each node u in a given set B t,k is then subject to a false-discovery rate correction as in Wilson et al. [2014] to obtain a set of significant nodes for each time point. However, our objective is to find sets B such that for each v ∈ B, v is significantly connected to u across all time points T , so we detail some additional steps to summarize and combine the sets of significant nodes. In the following sections 3.4.1 -3.4.3 we elaborate on the procedure to select significant bordering nodes across time. We also describe the initialization step at iteration k = 1 later in 3.5.1.

False Discovery Rate Correction
Given a test statistic S(u, B k,t , Gt), we derive time-ordered set of sets of p-values associated with nodes in graphs Gt obtained from (9), we use a procedure based on that of Benjamini-Hochberg [Benjamini and Hochberg, 1995]. For a fixed time t and fixed iteration step k, we fix set B k,t and find all the nodes that are significantly connected to it across all time t = 1, ..., T , and denote the output set as M k (B k ), which is defined in more detail in (12). Such a procedure is described as follows: Benjamini-Hochberg False Discovery Rate Procedure 1. Given a candidate set B k,t , compute p-values pu := p(u, B k,t ) for all u ∈ 1, ..., n in graph Gt as in (9). Note that nodes that do not border B k,t will have a p-value of 1.
2. Order n vertices of Gt such that the ordered p-values pu have ranks j(u) 3. Construct an adjusted p-value where at is the decay term contingent on the hypothesized trend described in the following Section 3.4.2 4. Define threshold τ (p) := max{pu : p * u ≤ α} , then obtain then the for all uj where Then update node-set mt(B t,k ) with the u's such that (11) holds. Following this step, we test for trend in Section 3.5, to be described later.

Time-Decay Adjusted False Discovery Rate Threshold
The Benjamini-Hochberg (BH) procedure uses the threshold cutoff of τ (p) to filter the sorted pvalues p * u . However, we augment the threshold of rejection contingent on whether the hypothesis of the interconnected communities are significantly increasing , decreasing, or stable in connectivity over time. We do so by introducing an exponential decay term to modulate for the shifting time-window of significance. Now we define the FDR adjuster at as : The terms a + 0 and a − 0 are defined such that at is 0 at time 1 and 1 at time T if the trend is increasing, and 1 at time 1 and 0 at time T if the trend is increasing: The time decay term within the exponential is normalized by T in our application in order to make the adjusters for the positive trends slightly convex and negative trends slightly concave with time in order to be slightly more permissive about selection of nodes that exhibit decreasing trends. Decreasing the adjuster slightly faster with time when the trend is hypothesized to be decreasing allows us to to be more aggressive about identifying clusters that drive down the overall ridership, to tailor to the hypothesis in Section 2 . A linear multiplier, or a different normalizer (such as n) to scale by the size of the graph), may also be appropriate in different settings.
If the trend amongst the connectivities of nodes is posited to be increasing, then the algorithm allows a more permissive selection of 'significantly' bordering nodes when time t is early. however, the algorithm is more penalizing (up to the rate of ordinary BH) when t approaches T . When time t is 1, then at is equal to zero, such that all p * u are zero and will automatically be counted as significant if u borders B k,t . When t is T , however, then at is 1 and so the adjusted p-value penalizes in the same way as in BH. The threshold for selection of the maximum allowable p-value increases as t decreases, so that nodes that only connect negligibly at earlier time points, but strongly at later points, are identified as significant.
Conversely, when the objective is identifying a decreasing trend, the exponential multiplier to the threshold is subtracted by one, so the same kind of correction is made in reverse. Because the multiplier to the adjusted p-value is always less than 1, the procedure is always less conservative than the Benjamini-Hochberg method and lets nodes that otherwise would not be significant at a given time-period be significant based on their potential to be significant given their trajectory. If the trend is posited to be neutral, then we use the ordinary BH rejection procedure.

Bonferroni Interval For Bordering Frequencies
To determine whether each node v is significantly connected to u at all times, we apply a second Bonferroni correction step. Such a correction is not applied to the p-values directly, but to the tabulation of frequencies of nodes whose p-values have been deemed significant by the previous correction. At a given iteration step k, For each B k,t across all time t = 1, ..., t, we use a product of Bonferroni confidence intervals to define significance of neighboring frequency . For arbitrary set B, we define mt(B) as the set of nodes that are found to be 'significantly bordering' described by Section 3.4.1: mt(B) = { number of vertices significantly bordering B at time t} we define F DI α,k to be the false discovery interval of all significantly adjacent nodes to node B k,t For each v ∈ B 0 k , we denote the bordering frequency Nv(B k ) as: Nv(B k ) = {number of times v is significantly bordering B k,t across time t}.
The neighboring frequency is significant if F DIT < Nv(B k,t ): that is, B k,t borders v enough times across t for it to be significant across all time. Finally, we deduce that v is significantly connected to B k across all time t and we take the union of all of these nodes and denote the set of these nodes as M k (B k ).
A significant Nv(B k ) suggests that v is more frequently bordered across time than other nodes that are identified to be 'significantly bordering' in 3.4.1. Large mt(B k,t ) for all t signifies a large collection of nodes that significantly border B k,t and results in a false discovery interval that is close to T , and hence v must border B k,t for nearly all time T for it to be significant. Conversely, if mt(B k,t ) is small, then the required frequency for v to be significant is not as high. During this step, we assume B k,t are independent under the null condition.
The resulting set represent the nodes that are significantly connected across time, given the appropriate time-window adjustments if the trend of clusters is postulated to be increasing or decreasing. Following this step of the algorithm, we check whether the trend of normalized connectivities is actually as hypothesized; this method is described in detail in the following section.

Significance Testing for Trends
Given an active set of nodes B k at iteration k and alternate between the calculation of significantly bordering nodes in (9) and significantly changing trends, which we describe in this section. At iteration step k, we first find all nodes v * that are significantly bordering across time as described in Section 3.4 and label the collection of these nodes as M k (B k ) as in (12). Following this step, we assess the significance of the addition of the candidate nodes v ∈ M k (B k ) to the existing time series of the active set B k by obtaining the test statistic for node-set connectivity S(u, B, Gt): We then find the time trend βv,B k using (14) for Before iterative testing in the general case for step k, we first initialize according to section 3.5.1.

Initializing Time Trend of ξ uv
To initialize the iterative search procedure, we examine all individual nodes u ∈ 1, ..., n. We calculate M0(u) for all B0(u) = u following the procedures from 3.4 at iterative step k = 0. Within M0(u), we calculate each normalized Wuv,t|Auv,t by the following equation for all v that are significantly connected to u across all time T : Next, we find the linear trends of each Zt(u, v) across time t = 1, ..., T and take the nodes with trends that are either significantly positive or negative. We write Z(u, v) as the time series of Zt(u, v) from t = 1, ..., T . The trend coefficient is calculated as the coefficient from a simple linear regression.
βuv is determined to be significantly increasing, decreasing, or stable (neutral) using the method described in the following section 3.5.2, but only using a single node v in place of a set B. If βuv is significant at the α level (measured by its t-statistic in a simple linear regression ), then we label the nodes v that are significantly connected and significantly increasing or decreasing with the initializing node u as v * * , and construct an initializing set B1 with these nodes {u, v * * } for the iteration step k = 1.

Energy Barrier for Equivalence Testing
Given that a time-series of networks yields significantly interconnected communities of nodes over time, we specify three different modes of relationships that compose the definition of intertemporal community. We posit that the three types of intertemporal communities of interest are those that are increasing in connectivity, those that are decreasing, and those that exhibit a neutral trend; these categories are specified before the application of the algorithm and as such every step of the procedure varies based on the clustering trend hypothesis in question.
To this end, we use the equivalence testing method [Schuirmann, 1987, Dixon andPechmann, 2008]. Even if a trend is significant under the assumptions of ordinary least-squares, it still may be regarded as negligible if its magnitude is small. We introduce a bounding parameter, or energy barrier U > 0, which serves as the threshold for a linear trend that is either positive, negative, or around zero. A positive U is chosen such that the significance of the trend is based on the difference between its estimate and U . Similarly, −U is used as a upper bound for the trend if it is negative. Finally, communities with a stable trend in connectivity are found by applying a test of trend equivalence with a symmetric bounding interval of [−U, U ] about zero. Pre-specifying a U before running the algorithm allows control over 'how large' the desired time-trends are in the the identification of intertemporal communities.

Increasing and Decreasing Trends
For the time trend expressed in (14) given a set B, node v, we test for hypotheses involving a symmetric interval [−U, U ] close to zero. These hypotheses approximates a null hypothesis of zero in equivalence testing, the null hypotheses are adjusted as follows: We calculate the significance of βv,B using the difference of the estimates as well as the (prespecified) upper and lower bounds of the trend. In order to filter out the trends that are negligible, we perform a t-test for the regression statistic subtracted by the upper or lower bound U , divided by the standard error of the estimate, s(βv,B). Defining such a bound allows us to exclude the very small but still significant trends and only find clusters that are sizeably increasing or decreasing.
The corresponding p-values of tupper and t lower , respectively, with significance α/2 (for one-sided tests) with degrees of freedom n − 2. Typical of ordinary least squares, the degrees of freedom are discounted by the slope and intercept terms.

Neutral Trend
To determine the t-statistic of a negligible trend, we utilize the approach outlined in Dixon and Pechmann [2008]. To test for whether a trend is negligible, the typical hypothesis test for a regression coefficient is inverted and split instead into two one-sided tests.
H0,a : βv,B ≥ U H1,a : βv,B < U, To test for such a hypothesis, Dixon and Pechmann [2008] used the following pair of t-statistics and obtained the corresponding p-values for the probability of the alternative hypothesis by taking the maximum of the p-values associated with the t-statistics t neutral,a and t neutral,b , respectively, with significance α and with degrees of freedom n − 2.

Iteration and Overlap Filtering Steps
After the procedures for selecting nodes that are both significant in connectivity (Section 3.4) and trend βv,B pending on the posited direction of trajectory (Section 3.5). We derive p-values from the t-statistic of the simple linear regressions corrected with the (regular) BH correction with significance level α/2 for the one-sided tests of increasing and decreasing trends, and α for the sum of the p-values associated with t neutral,a and t neutral,b for the test of neutral trend. We update the set B k+1 with the inclusion of the new nodes v that are both significantly connected to B k across all time t and has a significant trend according to the trend hypothesis. The iterative testing procedure is repeated until the set becomes stable, where B k = B k+1 for all candidate sets. Typically, this process takes 4 or 5 iterations.
After the iteration steps produce stable sets , the sets are filtered by their Jaccard overlaps in the same way as in Palowitch et al. [2018], He et al. [2019]. We use an overlap threshold of 0.50. After filtering by Jaccard overlaps, communities of size 3 or less are removed, as dyadic or triadic relationships between nodes may be too localized to be meaningful in a larger scale.

Results
We ran the algorithm under a range of values for α and U . Results under range of tuning parameters suggest similar geographical patterns of groups of bike stations that are increasing or decreasing in connectivity. Increasing clusters are mostly found in the Northeastern part of the city. Decreasing clusters are located in a pattern of a ring around the Downtown Loop, but most significantly in the Southwest, suggesting that decrease in overall ridership between 2016 -2018 is mostly accounted for by these clusters.

Primary Findings
We first discuss findings after fixing α at 0.05 and U ≈ .68. We fixed U at around .68 because it captured clusters of large sizes across all trend categories and also showed showed distinct geographical divisions. Under these tuning parameters, we found 5 clusters with decreasing connectivities over time and 5 clusters with increasing connectivities. Only 1 cluster was found to a neutral, or stable, trend at the 0.05 significance level. We call these communities of nodes decreasing, increasing, and neutral clusters for brevity.
There is a stark division between trends in ridership among the North and South parts of the city (see figure 4). When normalized by overall graph connectivity, and at the 5% significance level, clusters with significantly decreasing trends were mostly found in the South and West sides, while clusters with significantly increasing trends were all nearly found in the northern and central parts of the city. Interestingly, the decreasing clusters form a pattern that maps to a nearly concentric outer ring around the central part of the city, while the increasing clusters stretch from the Loop northwards along the lake shore.
The largest cluster is decreasing in connectivity over time and consists of thirteen stations, eight of which are in the Near West Side community area are five of which are in the Pilsen neighborhood (in the Lower West Side community area). The Near West Side is a large community area that is mostly affluent but also containing many diverse neighborhoods, such as Little Italy and the University of Illinois in Chicago. The Lower West Side is a much less affluent and predominantly Latino neighborhood, although its demographic composition is rapidly changing [Betancur and Kim, 2016].
The neighborhood that houses the most stations that are in increasing clusters is in the historically affluent Lincoln Park. Most of the clusters contain stations in Lincoln Park. Interestingly, three stations in Lincoln Park connect to the relatively distant the Western and Monroe station in the Near West Side, which is the only station on the west side that is an increasing cluster. There is only one cluster located in the central part of the city that has a neutral trend. Three of the stations in this cluster are located in the Loop, which is the central business district. Such stability may be due to consistently high volumes of usage in the business areas of the city.

Sensitivity to Parameters
The sensitivity of cluster detection to the choice of parameter is often overlooked [Austwick et al., 2013]. The significance parameter α and the bounding parameter U produce generally similar results over a range of values ( Figure 5). Shifting U to 0.005 while keeping α fixed at 0.05 induces discovery of more increasing and decreasing clusters, but the bound is too tight for any significant sets of to found under the hypothesis tests in (16). The clusters found under these settings generally agree with the geographical trend in 4. Negative clusters make an even more apparent concentric circle around the central part of the city, but more clusters are now found in the northwestern neighborhood Logan Square as well as in Woodlawn and Hyde Park in the South East Side. More increasing clusters in the central and northeast part of the city are found when U = 0.005. A cluster comprising of stations such as Adler Planetarium and Millenium Park in the Downtown Loop, which connects many cultural and tourist amenities was found under this barrier. Several other clusters in the Near North Side and the Loop were found under this threshold that were not under the settings in 4.
In general, we see a trade-off between finding increasing or decreasing clusters and stable clusters hinging on the choice of U . If U is larger, then there is more room for a trend to be classified as stable, but less so for increasing or decreasing trends. The discovery of more increasing and decreasing clusters when U is set at 0.05 compared to U = 0.07 suggests that these clusters are not increasing or decreasing in connectivity at the same rates as the other clusters. When U is further raised to 0.01, while α is still fixed at 0.05, only few clusters remain in the increasing or decreasing categories, but more clusters whose connectivities are stable over time are present.
When U = 0.01 (row 2 of figure 5), the three decreasing clusters are in Near West Side and Lower Side, Armour Square, and Bridgeport.The clusters that are significantly increasing are in the Near North Side, Lakeview, and the Loop. These clusters are increasing and decreasing at more larger consistent rates than others that were found with lower barrier. However, because α is the same as in row 1 of 5 , they are subject to the same connectivity.
The neutral clusters are mostly also in the Loop and the Near North Side. These clusters share a lot of commonalities with increasing clusters at less stringent thresholds, so it is reasonable to assume that they probably have a moderate increasing trend. The clusters mostly are located in regions with high foot traffic around downtown, museum district, and lakefront.
Row 3 of Figure 5 shows one cluster that was found when U is set at 0.001 and α is decreased to 0.01. All stations of cluster are found in the neighborhood of Bridgeport. This cluster has a lower barrier than results in Figures 4 and 5, but has a higher significance threshold, so these nodes are more interconnected in later time-points compared to other decreasing clusters, but has a slower rate of decrease. Figure 6 shows more in-depth behavior of these 4 stations in Bridgeport. The top left figure Sum of Trips shows that the total sum of edges between all four of the of the members of the community starkly decreases across time. Taking the 'raw' sum of the trips between the four stations shows very strong seasonal signals. However, by normalizing each time unit by the nodes' degrees and strengths using the procedure detailed in 3.4, the total normalized connectivity sum Z(B) exhibits a steadily decreasing trend, as do all of the normalized node-set scores Zt(u, B).
It is perhaps important to note that the interaction between α and U isn't entirely linear or monotonic. Though a decrease in α intuitively corresponds to an increase in U , and empirically they may be related in this way, it is not always the case that . We can see in figures 4 and 5, that as U becomes larger and α stays the same, different clusters that didn't exist at a lower threshold may appear. This is because of the Bonferroni and BH steps, a lower barrier U may cause more nodes to be found significant at the BH reject, but then results in a lower trend with the inclusion of these nodes.

Conclusions
The method devised in this study can be applied to a variety of data in network time-series format, particularly human mobility networks. In depth analysis of neighborhood changes may be an immediate extension of the substantial aspect of this study in gleaning demographic trends. In particular, changes in neighborhood ethnic compositions across years may be impactful. We only use a snapshot of demographic data by community in 2015. The method can be applied to bikeshare networks in other cities to compare the relative rates of growth and decline in bikeshare usage. We will elaborate on the theoretical properties of intertemporal community detection in future work.
The proposed method yields several substantial and methodological contributions. We applied a novel community detection method on time-series of networks and found clusters that, when validated with visual inspection as well as neighborhood demographics, yield meaningful results. Increasing clusters mostly are contained within affluent, mostly white neighborhoods. Neutral trends tend to be around the downtown loop area, while decreasing trends form a concentric pattern along the edge of the downtown area. These geographical patterns suggest that there may be latent factors undergirding the decreasing trends that are driving the overall city-wide decrease in usage.
There may be several explanations for this underlying signal that cause the clusters to decreasing connectivity. One reason may be that because the neighborhoods housing the stations that form decreasing clusters are slightly less affluent than those that are increasing, hence leading to diminished usage. Another reason may be that decreasing trends may be symptoms of displacement in destabilizing long-term inhabitants and in turn, steady ridership. The fastest-gentrifying neighborhoods that are mostly non-white are also subject to decreasing trends in connectivity. Neighborhoods that are being gentrified are subject to higher rates of displacement, and in turn, loss of activity. Our study may be descriptive of such processes.
If we assume that people use Divvy at similar rates across demographic stratifications and that most of bikeshare usage occurs around residents' immediate neighborhoods, then bikeshare usage parallels general 'outside activity' and human mobility across neighborhoods. Though this assumption maybe generalizing, it is one framework which justifies the possible interpretations of our findings in light of justifying the connectivity trends between stations around the edge of the city center.
Most decreasing clusters which have stations in gentrifying neighborhoods are also connected to a more affluent neighborhood. These relationships are seen in the clusters that straddle Pilsen and the Near West Side, Bridgeport and Armour Square, and Wicker Park and Logan Square. Figure 4 shows a stark geographical divide between decreasing and increasing clusters, in which decreasing clusters form nearly a circular outer ring around the city, but increasing clusters are all concentrated in the North East side along the lake, which do tend to map to white and historically affluent neighborhoods.
Most decreasing clusters are in the fastest-gentrifying neighborhoods. That the clusters are decreasing may be due to construction projects in neighborhoods exhibiting economic growth; con- struction projects make it more difficult for people to ride bicycles thus less likely for people to use bikeshares. Additionally (or alternatively), decrease in connectivity may be indicative of displacement in gentrified neighborhoods. In particular, the decreasing clusters may show the effects of greater volumes of displacement amongst existing residents than emigration of new residents. Our results show that gentrifying and recently gentrified neighborhoods are experiencing a steady rate of (relative) decline in Divvy ridership. Prior research has shown that when areas become gentrified, even as new residents move in, inhabitants who are pushed out leave in much greater numbers [Betancur and Kim, 2016] and new housing typically becomes less dense [UC Berkeley Health Impact Group , 2010]. Moreover, even a short dip in vacancy can disrupt long-term usage of bikeshares. If a resident who consistently uses bikeshares in their neighborhoood moves out and someone else immediately moves in, they may not adjust to the routines that the previous inhabitant has established.
We proposed a novel method to cluster human mobility networks that vary across time. Our method combines usage of a null model in relation to network topology as well as a null model in relation to the expected trajectory of the graph evolutions. By adhering to a configuration null model for every time point, we are able to remove some of the global (city-wide) seasonal effects and focus on mapping the latent relative connectivities between Divvy stations. The methods used in this study may be applied to other situations where it is important to study the evolution of structures within networks.

Availability of Data and Material
All the data are available on the Divvy website [Divvy, 2019] (https://www.divvybikes.com/system-data). The code to implement the methods described in this manuscript is available upon request.