Computer Science Mathematics
Computer Science-2021.gz:143353
Computer Science
b34454dfa563a05f4fed792dd539ad38062d1d5f
Geographical window based structural similarity index for origin-destination matrices comparison

Abstract Most traditional metrics compare origin-destination (OD) matrices based on the deviations of individual OD flows and often neglect OD matrix structural information within their formulations. Limited metrics exist in literature for the structural comparison of OD matrices. One such metric is mean structural similarity index (MSSIM) that computes statistics on groups of OD pairs defined by local sliding windows. However, MSSIM can result in different values based on the choice of the size of the window. In literature, no clear consensus has been reported on the level of acceptability of the window size and the resulting MSSIM values. Addressing this need, we propose the concept of geographical window, and develop geographical window based structural similarity index (GSSI) that exploits OD matrix structure by computing statistics on the group of OD pairs that are geographically correlated. Compared to traditional sliding window based MSSIM, the advantages of GSSI technique identified from real case study application are (a) it preserves geographical integrity; (b) it compares results with physical significance; (c) it captures local travel patterns; (d) it compares large-scale sparse OD matrices; and (e) it is computationally efficient. A thorough sensitivity analyses suggest that GSSI is a robust statistical metric and has potential for practical applications such as, benchmarking different OD estimation methods; improving the quality of solution by maintaining structural consistency in the OD estimation process; and identifying gaps in the transit service by comparing local (within a geographical window) travel patterns of car and public transit.


Introduction
An origin-destination (OD) demand matrix represents the distribution of travel demand between different OD pairs within a transport network. It is a key input into most traffic simulation models and it essentially provides two major aspects of travel information: first, the cell values of the OD matrix represent flows between individual OD pairs (we can also term it as the mass); and second, the arrangement of OD pairs reflects the underlying structure (or skeleton) of OD matrix.
In the OD estimation process, knowledge of OD structure helps to maintain structural consistency in the OD estimates (Behara et al., 2020b;Doblas & Benitez, 2005). Similarly, accounting for the OD structure in the comparison of the OD matrices from different solution algorithms (on a benchmark network) provides a holistic quality measurement (Djukic et al., 2013).
In literature, different methods (details below) are used to quantify the similarity between two OD matrices. If the structure of the OD matrices is also considered in the similarity estimation, then we term it as structural similarity. Two OD matrices have perfect structural similarity if their structures are similar with zero differences in the OD flows.
Some notable traditional measures widely used to quantify similarity of OD matrices are: root mean square error (RMSE) (Frederix et al., 2014;Hellinga & Aerde, 1994); normalized root mean square error (RMSN) (Antoniou et al., 2004;Jaume Barcel o et al., 2013;Frederix et al., 2014); mean square error (MSE) (Cascetta, 1984); mean absolute error ratio (MAER) (Kim et al., 2005); mean absolute percent error (MAPE) (Nigro et al., 2018); goodness of Theil's fit (GU) (J. Jaume Barcel o et al., 2013); R-squared (R 2 ) (Tavassoli et al., 2016); and entropy measure (E) (Ros-Roca et al., 2018). Refer to Hollander and Liu (2008) for more information about the traditional metrics. Although the formulations of these metrics are mathematically simple they often fail to capture the structural differences between OD matrices (Djukic et al., 2013). Thus they are incapable of accounting the OD structural properties such as trip productions/attractions, destination choices, geographical correlations (Antoniou et al., 2016). Refer to example in the next section that demonstrates the inability of traditional metrics to capture OD matrix structural differences. Also, these standard measures are incapable of recognizing the spatial and temporal proximity of flows in the adjacent cells . Total demand scale (TDS), proposed by Bierlaire (2002), is another measure to assess the quality of estimated OD matrices. It measures intrinsic under-determinacy of OD estimation problem by computing the difference between minimum and maximum levels of total demand of the corresponding OD estimates. Lesser the difference between the minimum and maximum OD combinations, better is the quality of estimated OD. However, TDS is not a suitable approach for applications other than OD estimation problem (e.g., comparison of OD matrices from different day types).
Limited metrics such as mean structural similarity index or simply MSSIM (Djukic et al., 2013), 4D-MSSIM (Van Vuren & Day-Pollard, 2015), Wasserstein metric (Ruiz de Villa et al., 2014), and normalized Levenshtein distance for OD matrices (NLOD) (Behara et al., 2020a) are available in the literature that could account for OD matrix structural information. These metrics are popular in other disciplines and are introduced in transport applications as discussed below.
The metric MSSIM computes statistics i.e., mean, standard deviation and covariance comparison on groups of OD pairs defined by local windows, and can identify structural differences better than traditional metrics such as R-squared and GEH (Day- . Refer to set of Equations (2) to (3) in the next Section for mathematical formulations. Originally, MSSIM is applied for structural comparison of natural images considering that neighborhood pixels are correlated. However, the major drawback is that MSSIM is sensitive to window size and can result in different values based on the choice of local window size (i.e., size of OD pairs' group). Additionally, application of MSSIM on OD matrices needs to account more attributes that are OD related. For instance, (a) the correlation between OD pairs can exist for a variety of reasons such as spatial proximity, and trip distance, and there is no guarantee that all OD pairs within the window are correlated; (b) the stability constants (refer to Equations (2) to (3)) used in its formulation could be network specific; and c) the weightages to SSIM components need to be further explored especially while comparing sparse and dense matrices (Day- . More details about the formulation are provided in the subsequent sections.
To define a group of OD pairs that are in spatial neighborhood, Van Vuren and Day-  proposed 4D-MSSIM. It is an extension to MSSIM that additionally accounts for spatial proximity between OD pairs through Euclidean distance. The nearby OD pairs need to be computed for every OD pair in turn. The proximity estimation is only a onetime computation for a particular study area, and the idea of 4D-MSSIM is more in tune with OD matrices. However, it was still in an explorative phase as Van Vuren and Day-  recommended further testing its implementation on a range of networks under wider set of circumstances.
Wasserstein distance is defined as the minimum total travel time required to assign the trips between OD pairs of matrix "X" using an assignment compatible with matrix "Y." It is a linear programing problem that yields distance in terms of travel time or minutes per driver. This metric is based on optimization technique and is computationally more intensive than MSSIM or 4D-MSSIM. Refer to Ruiz de Villa et al. (2014) for more details. Similar to Wasserstein, normalized Levenshtein distance for OD matrices (or NLOD) also compares OD flows through an optimization approach. It defines OD structure as the preference of destinations from each origin. This approach is inspired from the traditional Levenshtein widely used in computational linguistics such as comparison of string sequence. Refer to Behara et al. (2020a) for more details.
Research on structural similarity of the OD matrix is limited. This paper aims to rationally extend the applicability of SSIM on OD matrix by considering the geographical boundaries. The motivation for which is twofold. First, there is lack of consensus on the size of the window for SSIM. Second, the transport analysis zone (TAZ) boundaries can be aggregated at different levels with similar socio-economic and land use pattern.
The increasing need to minimize the gap between land use and transport modeling motivates the modelers to standardize the statistical boundaries and TAZ boundaries. For instance, Australian Bureau of Statistics (ASGS, 2017) defines hierarchy of geographical areas for the release of statistical information. This includes statistical areas (SA) for four levels: Statistical area level 1 (SA1) to statistical area level 4 (SA4). Generally, SA1 has a population of between 200 and 800 persons; SA2 level normally reflects the sub-urban level and are aggregation of SA1s; SA3 is designed at the regional level and is an aggregation of SA2s; and SA4 reflects labor market within each State and Territory and is an aggregation of SA3s. For transport modeling the SA boundaries are being considered as TAZ with certain modifications where needed (Weston, 2019).
We propose to compute SSIM on the geographical windows that are defined by geographical boundaries and term this extension as geographical window based structural similarity index (GSSI). The sensitivity of proposed index is thoroughly tested with a case study application using Bluetooth based OD (B-OD) matrices from the Brisbane City Council (BCC) region, Australia.
The remainder of the paper is structured as follows: First, the study demonstrates the need for OD matrix structural similarity; and then presents the traditional MSSIM method. Next, it explains the development of the proposed measure -GSSI followed by discussion on the advantages of GSSI over MSSIM with a case study application on Brisbane city. The study then tests the robustness of proposed method through sensitivity analysis; discusses about its performance; and finally, concludes with future research directions. The Appendix section at the end of this paper provides a conceptual comparison of GSSI with other structural (dis)similarity metrics.

Need for OD matrix structural similarity
The OD flows between different OD pairs can have strong correlations due to sharing of similar structural properties such as trip production/attraction, trip purpose, spatial proximity, travel cost, and similar routes. While the temporal correlations exist among flows for the same OD pair (Nigro et al., 2018), the spatial correlations refer to the relationship between OD flows from different OD pairs during the same time period (Djukic, 2014).
Consider an example shown in Figure 1 for the demonstration of spatial correlation between OD pairs. In this example, we consider trips between OD pair/s represented at different zonal level hierarchy; that is, SA-3 and SA-2 for higher and lower levels as shown by thick and thin grey colored boundaries in Figure 1. The SA3 origin "O" and destination "D" are Nathan and Brisbane Inner (as highlighted by yellow colored boundaries), respectively connected through Pacific Motorway. Similarly, the trip boundaries at SA2 level are origins (o 1 , o 2 , o 3 , and o 4 correspond to Tarragindi, Moorooka, Salisbury-Nathan, Coopers Plains, and Robertson, respectively); and destinations (d 1 , d 2 , d 3 , and d 4 are Brisbane City, Spring Hill, Fortitude Valley, and New farm, respectively).
In Figure 1 we can see that SA2 OD pairs (i.e., at lower zonal level) could be correlated when they share similar structural properties. For instance, the SA2 zones (say, o 1 to o 5 ) are geographically correlated because they belong to the same SA3 zone (i.e., Nathan); they might have similar trip generation characteristics; and mostly share the same route (i.e., Pacific Motorway) to reach destination "D." The traditional metrics compare OD matrices based on deviations between individual OD flows and thus cannot recognize the spatial or temporal proximity of movements in adjacent cells and the complex OD structural information . To demonstrate this, let's consider an example of comparing OD matrices M i (5M 1 or M 2 ) with a reference OD matrix M R (Figure 2). Here, M 1 is simply 1.1 times M R , and M 2 is chosen randomly. Traditional metrics namely MSE, RMSE, GU, and MAE are used for matrix comparison as shown in the Equations (1)-(1c) where M R, w and M i, w are the reference and query OD flows of w th OD pair from OD matrices M R and M i , respectively (total number of OD pairs in both matrices is W ¼ 16 each). The results of comparing matrices M 1 and M 2 with M R are presented in Table 1. The first column of Table 1 presents the metrics. The second and third column are the values from metrics for both cases, respectively. M 1 is structurally similar to the reference matrix than that of M 2 . The traditional metrics fail to quantify the structural differences as the values of the  metrics are approximately same for the two comparisons.
From the above example, we could see that traditional metrics are unable to account for the OD matrix structural properties. This demands further need to focus research on structural (dis)similarity measures as additional metrics in transport applications that involve OD matrices comparison. One such metric is mean structural similarity index (MSSIM) which could potentially differentiate the above two matrix pairs. The application of MSSIM formulation (discussed in detail in the next Section) on the matrices illustrated in Figure 2  Mean structural similarity index (MSSIM) Wang et al. (2004) proposed mean structural similarity index (MSSIM) as a quantitative measure to compare the quality of two images. They highlighted the inability of traditional metrics in accounting the (visual) structural differences between the images. In their paper (refer to Figure 5 in Wang et al. (2004)) they showed that two images had the same MSE but different SSIM when compared with the same reference. Djukic et al. (2013) argued that SSIM can be applied to compare two OD matrices because OD pairs are analogous to pixels of an image (refer to Figure 6.2 in Djukic (2014)). The details of the SSIM formulation are presented in the following sub sections.

Local sliding window and MSSIM formulation
The traditional SSIM computes statistics on the local window (consisting group of pixels or OD pairs) that slides pixel by pixel (cell by cell) over the entire image (OD matrix). It is often referred as "local sliding window" and is generally of size m Â n (where, m Â n ( M Â N i.e., size of OD matrix). The concept of sliding window is used to enhance the ability of SSIM to account image distortions by computing local statistics (Brooks et al., 2008). For ease of explanation consider an example as presented in Figure 3. Here, two 4 Â 4 OD matrices, X and Y, are presented in column one and two, respectively. These two OD matrices are  compared using SSIM. The local sliding window of 2 Â 2 sub-matrix (m ¼ n ¼2) is considered and is represented as colored cells. This window slides cell by cell over the entire OD matrix and in the current example we have 9 pairs of sub-matrices (refer Figure  3a-l). The local SSIM computes the structural similarity between the sub-matrices corresponding to the windows from both OD matrices. The MSSIM is computed by averaging the local SSIM values computed for all sliding windows. The formulation for local SSIM is the product of three individual formulations (Equations (2)-(2b)) related to mean, standard deviations and coefficient of correlation between the groups of OD pairs. where X and Y represent two OD matrices to be compared; x and y represent the group of OD pairs within local windows in both matrices. The individual components of L x s , y s ð Þ and C x s , y s ð Þ compare the mean values (l x s and l y s ) and the standard deviations (r x s and r y s Þ among OD pairs within the same (s th ) local windows from both matrices. The component STR x s , y s ð Þ compares the structure by computing coefficient of correlation between s th group of OD pairs in both matrices. The constants C 1 , C 2 and C 3 are meant to stabilize the result when either mean or standard deviation is close to zero. Generally, C 3 is assumed to be C 2 =2: Previous studies suggested values of 10 À10 and 10 À2 for C 1 and C 2 , respectively (Pollard et al., 2013). In the analysis conducted for this study, the OD flows within the SSIM window are not all zero, hence we assume both C 1 and C 2 as zero. The parameters a, b and Y are used to adjust the relative importance of mean, standard deviation, and structural components, respectively, and are generally assumed to be equal to 1. The SSIM x s , y s ð Þ (Equation (3)) is the structural similarity of s th local window (of dimensions m Â n) from both matrices; and mean SSIM or MSSIM X, Y ð Þ (Equation (3a)) reports the overall structural similarity of the OD matrices, X and Y that is computed by taking average of SSIM values from S sliding windows. The values of both SSIM and MSSIM lie between À1 and 1. The value of 1 implies that matrices are the same while the reverse is true when value is À1. In Figure 3a, the SSIM value for local window is 0.5963 and MSSIM computed over 9 local windows is 0.6777.

SSIM extension in literature: 4D-MSSIM
We know that consecutive zonal numbering (i.e., 101,102, 103 and 104 in Figure 3) in OD matrices does not always guarantee geographical adjacency between the zones. If the zones are not adjacent, then spatial adjacency among OD pairs is further not guaranteed. To apply the SSIM's rationale of "near-by pixels" in the OD matrix context, spatial proximity between OD pairs need to be considered. In this direction, Van Vuren and Day-Pollard (2015) developed 4D-MSSIM that accounts for spatial proximity among nearby OD pairs. The proximity is computed using Euclidean distance between the spatial coordinates of OD pairs. The Euclidian distance (d s b a , s d c À Á ) from OD pair s b a (flow from zone a to zone b) to OD pair s d c (flow from zone c to zone d) is expressed using Equation (4) where the coordinates of the centroid of zones a, b, c, and d are (x a , y a ), (x b , y b ), (x c , y c ), and (x d , y d ), respectively.
Contributions of nearby OD pairs are calculated using Gaussian weights (c) expressed as shown in Equation (5). Here, "x" refers to Euclidian distance between OD pairs and "r" is the network attribute that is generally average trip length or diameter of sample of zones in the study area. The shape of Gaussian plot with weight (c) in y-axis and Euclidian distance in x-axis is shown in Figure 4.
The final 4D-MSSIM is the MSSIM computed on the group of weighted OD flows. While the concept is interesting, it was still in an exploratory phase and there were few questions that needed to be answered. First, the cluster of "near-by" OD pairs is sensitive to the selection of initial OD pair; and second, Euclidian distance might not be an accurate measure of spatial proximity between OD pairs. For example, travel distance between two sets of OD pairs separated by a natural barrier (say, river) could be much higher than a mere Euclidian distance. In regards to Equation (5), Van Vuren and Day-Pollard (2015) mentioned that "The value in the denominator deserves some attention; it would make sense to link this to attributes of the network e.g., the average trip length or average zone diameter and would thus be mode and possibly purpose-specific. We have not yet analyzed the effect of different assumptions in this respect."

Need for improvements in the SSIM approach
Although the formulation of SSIM has the ability to capture the structure of OD flows, the approach adopted needs further improvement due to the following reasons: First, SSIM is sensitive to the size of local window due to which no clear consensus has been reported on the level of acceptability of sliding window size and the resulting MSSIM values. To demonstrate the sensitivity of SSIM toward window size, consider MSSIM values computed using different window sizes (from 3 Â 3 to 20Â 20) for Monday and Sunday, and Monday and Tuesday Bluetooth based OD (B-OD) matrix pairs that are constructed from the BCC Bluetooth data (more details about B-OD matrix are presented in the next section). Figure 5 presents the results where blue line is for Monday and Sunday, and orange line is for Monday and Tuesday comparison. The x-axis represents the size of the local window and y-axis shows the MSSIM value. The order of OD pairs is same in the matrices for Sunday, Monday, and Tuesday. It is observed that as the size of sliding window increases, sensitivity of SSIM toward subtle structural differences between the OD matrices decreases. The MSSIM values increase as the sliding window size increases. Similar results are observed by Brooks et al. (2008) in comparing images using different window sizes. The rate of increment of MSSIM values is less for Monday and Tuesday pair as compared to Monday and Sunday pair. This is due to similar travel patterns between Monday and Tuesday (both working days) and less similar patterns between Monday and Sunday pair.
To circumvent this ambiguity about local window size literature suggests that SSIM should be computed over the entire OD matrix (Djukic, 2014). However, by doing so it will result in a statistical estimation that is less sensitive to structural changes within the OD matrix. According to the law of large numbers, the variance within a sample tends to decrease if the sample size increases. Larger window dimensions imply a greater number of OD pairs to be compared. Thus, the variance (distortion), and covariance (correlation distortion) parameters that capture structural changes within, and between OD matrices shall be reduced. Therefore, larger window size makes SSIM to be less sensitive to structural distortions, and thus, a smaller sized window is recommended.
Second, as the size of sliding window become smaller, the computational cost also increases. For instance, the number of pairs of local sliding windows that need to be compared from both OD matrices are (M-mþ1) Ã (N-nþ1).
Third, not all OD pairs within a local window are correlated, and the local SSIM value does not have any physical meaning attached to it. For example, the local SSIM value of 0.5963 in Figure 3a does not convey any physical significance if the OD pairs are not : o(i) : : : : : : : : : Geographical window representing group of lower zonal level OD cells from higher origin zone O(1) to higher destination zone D(2) Figure 6. An example to illustrate the proposed geographical window-based approach. correlated within the local window. Djukic (2014) tried to address correlations among the OD pairs from their flow values (especially if volumes are high) by matrix-reordering (i.e., sorting each row of OD matrix in the order of OD pair volumes). Thus, correlated OD pairs can lie in the same neighborhood i.e., all high-volume OD pairs on one side and remaining on the other side (similar to the arrangement of pixels in any natural image). However, comparison of rearranged matrices might not be appropriate if the order of both matrices is different; for example, comparing Monday and Sunday OD matrices.
To this end, the study proposes an extension to MSSIM through the concept of geographical window as further discussed in the following section.

Proposed geographical window based structural similarity index (GSSI)
In the proposed approach, we first arrange the origins and destinations of OD matrix in order of the geographical identity, and subsequently define the windows for SSIM analysis that are consistent with the geographical boundaries. Here, the window size varies with the geographical boundaries considered in the rearranged OD matrix. This is different from the traditional SSIM application where the size of the window is fixed. The local window with a defined geographical boundary is termed as geographical window.
Let us explain the concept of geographical windows with the help of Figure 6 where each cell of the OD matrix represents a lower level OD pair (o i to d j ). Here, the OD matrix is rearranged so that the lower level origins (rows) and destinations (columns) can be grouped into respective higher-level origin (O k ) and destination (D l ) zones. The number of higher zonal origins and destinations are K and L, respectively. For instance, in Figure 6, o 1 to o i from O 1 level is arranged together. The higher-level boundaries now define the geographical windows. The orange-shaded region represents a geographical window covering OD pairs from O 1 to D 2.
The STR and SSIM are computed over the geographical windows using Equations (2b) and (3), respectively. Similarly, the mean of STR and SSIM values computed over G geographical windows are referred as GSTR (Equation (6)) and GSSI (Equation (7)), respectively.
where X and Y represent two OD matrices to be compared; x g and y g represent the group of OD pairs within g th local geographical window in both matrices.
Other parameters are already explained in the earlier section. The values of both GSTR and GSSI lie between À1 and 1. Note that GSTR accounts only the structural component and GSSI compares overall structural similarity of OD matrices.

Demonstration of geographical windows with an example from Brisbane city
The study region that is considered for demonstration purpose is the BCC region (refer Figure 7). Several studies in the past have used Bluetooth in transport applications (Bhaskar & Chung, 2013;Chen & Bierlaire, 2015;Michau, Pustelnik, et al., 2017;Yang & Wu, 2018). The current study also uses Bluetooth travel information for the demonstration purpose. Figure 7 shows 845 Bluetooth MAC scanners (BMS) overlapped on SA3 zones of the BCC region. The trips observed as a sequence of BMSs are used to develop OD flows between scanner pairs (refer to ) for more details about deriving Bluetooth-inferred trips). The OD flows at BMS level are further aggregated to SA3 level based on the concordance between BMS and SA3 zones. Thus, we have SA3 level B-OD matrices of 20 Â 20 dimensions. Since they are aggregated trips, the B-OD matrices are dense (as shown for Monday in Figure 8a). It is to be noted that the Bluetooth trip ends do not represent actual trip ends of the Bluetooth equipped vehicles, and B-OD matrices are used only for the demonstration of proposed methodology. Figure 8 demonstrates the application of SA4-based geographical windows for comparing SA3 (20 Â 20) B-  (Figure 8b). B-OD matrices are chosen randomly and belong to a typical week. The 7 th March (Monday) was a typical working weekday and 13 th March (Sunday) was not a part of any long weekend. These two days in the week are chosen because it is well known that weekend travel patterns are different from that of a weekday, and ideally, the comparison between the two should show dissimilarity. The SA4 zones used in designing geographical windows are Brisbane East, Brisbane North, Brisbane South, Brisbane West, and Brisbane Inner. For example, consider a geographical window of SA4 OD pair Brisbane East and Brisbane North. It consists of SA3 OD pairs e.g., 30101 to 30201, 30202, 30203, and 30204; 30103 to 30201, 30202, 30203, and 30204. These SA3 OD pairs are geographically correlated because they belong to the same SA4 origin (Brisbane East) and SA4 destination (Brisbane North). The size of local geographical window corresponding to this SA4 OD pair is 2 Â 4 because Brisbane East and Brisbane North consist of 2 and 4 lower level (SA3) zones, respectively.
In the above example, the total number of geographical windows considered is equal to the number of higher order OD pairs, which is 5 Â 5 ¼ 25, and the GSSI for Monday-Sunday matrix comparison is 0.7231.

Advantages of GSSI over MSSIM: Brisbane case study examples
The advantages of the GSSI approach over traditional sliding window based MSSIM are discussed below in further detail with case study examples from the BCC region.

Structural comparison of local travel patterns
The concept of geographical window ensures geographical integrity and captures spatial correlation by computing statistics on all local (lower zonal level) OD pairs belonging to the same higher zonal level OD pair. Thus, it provides opportunities to compare the local travel demand distribution (travel patterns) between different suburbs of a region that a sliding local window is not capable of. From Figure 9, it can be illustrated that Sunday travel patterns differ majorly for suburb pairÀBrisbane South to Brisbane North. This is reflected from local SSIM value of 0.4653 (see Figure 9 (left) and bold value in Table 2). On the other hand, for another suburb pairÀBrisbane South to Brisbane WestÀSunday travel patterns are similar (if not exact) to that of Monday with a value of 0.8037 (see Figure 9 (right) and bold value in Table 2).
Analysis of local travel patterns during different day types can help policy makers to identify suburbs for better transport planning. For example, assuming that B-OD (car OD) matrices in Figure 8 are reliable, the travel demand distributed between Brisbane East and Brisbane West during Monday and Sunday are very similar with SSIM equal to 0.9517. Comparison of such local demand distributions with public transit patterns can help transport planners to identify gaps in the transit service. For instance, consider the suburb pair Brisbane East and Brisbane West. If SSIM (Monday, Sunday) for transit OD is far lower than 0.9517 (of car OD), then it implies that the transit service between the suburbs might need further improvement.

Comparison of window sizes
The size and shape of a geographical window is defined by the number of lower zonal OD pairs present within the respective higher zonal OD pair. Therefore, in this approach the local geographical window need not always be a square matrix and might not be of fixed dimensions. For instance, in Figure 8, the geographical window are of different sizes such as 2 Â 4 for Brisbane East and Brisbane North, and 6 Â 4 for Brisbane South and Brisbane West.
Although the dimensions of the window may vary, the GSSI values are found to capture structural distortions as effective as the small window size based SSIM does. To demonstrate this, we have compared a typical working Monday OD matrix (i.e., average of 45 working Monday OD matrices from 2016) with 40 B-OD matrices each from Saturdays (Figure 10(a)) and Sundays (Figure 10(b)), and with nearly 45 OD matrices from rest of the weekdays i.e., Tuesday, Wednesday, Thursday and Friday (refer Figure  10(c-f)). It is observed that GSSI plot closely align (same variation) with that of a 2 Â 2 sliding window for weekends, and to that of a 3 Â 3 sliding window for weekdays. Figure 10 demonstrate 12 different plots each (11 of which correspond to sliding windows of sizes ranging from 2 Â 2 to 20 Â 20 (dotted lines) and geographical window (continuous line)). For each plot, x-axis corresponds to different daily OD matrices and y-axis reflects MSSIM/GSSI values.

Ability to capture subtle structural distortions
Here, we compare GSSI and MSSIM with a real case study example. For the study area shown in Figure 7, we consider B-OD matrices from the days 7 th -13 th Mar (entire week), 10 th Aug (Queensland Ekka festival), 25 th Dec (Christmas Sunday), and 27 th Mar (Easter Sunday) in 2016, and compare each daily OD matrix with the reference Monday OD matrix (i.e.,  Mar 7, 2016) using GSSI/MSSIM and GSTR/STR as shown in Figure 11 and Figure 12, respectively. Note that in this application, 25 geographical windows (as shown in Figure 8) are used for GSSI/GSTR, and only a single window (of size equal to B-OD matrix) is used for SSIM/STR because the order of OD matrices are different for different day types. The similarities in results of GSSI and SSIM comparisons are as follows: 1. The OD matrices during the weekdays are structurally more similar to the Monday OD matrix than rest of the day types. This is reflected in higher GSSI/SSIM (> ¼0.9899/> ¼0.9950) and GSTR/STR (> ¼0.9938/> ¼0.9995) values for the weekdays in Figure 11. A closer observation of GSSI/SSIM values in Figure 11 revealed that travel patterns during Tuesday and Wednesday (red rectangle) were slightly different from Thursday  and Friday (blue rectangle). This shows that GSSI/SSIM are sensitive toward structural differences within the working weekdays. 2. Both GSSI and SSIM identified that travel patterns during weekends and rest of public holidays are very different from that of weekdays. Among the non-weekdays, GSSI/SSIM values for Saturday-Monday showed higher proximity because the percentage of work-related trips occurring during Saturdays is higher as compared to Sundays and rest of public holidays (Bhat & Misra, 1999).
Despite the above similarities, we have found that SSIM is less sensitive to the structural distortions among the non-weekdays. This can be concluded from the section of plots corresponding to non-weekdays, especially Ekka, Saturday and Sunday, in Figure   13. We can observe that GSTR values (green oval) could differentiate more than STR does (orange oval) for these day types. To understand these differences, comparison of local travel patterns between different sections of the network is necessary which is only possible with the GSSI technique. Refer to Figure 13 that helps to visualize these local structural differences among these day types. The rows and column headers in Figure 13 are SA3 origins and destinations. The Ekka festival generally happens at Brisbane show grounds in Bowen Hills (Brisbane Inner). Due to this, the trips attracted to it and its neighboring suburbs such as Newstead-Bowen Hills, and Newmarket in the North; Ashgrove, and Milton in the West; and CBD, Valley, and South Bank in the Inner city have different demand distributions on the Ekka public holiday as compared to usual Saturday and Sunday.  Figure 11. Comparison among daily OD matrices using GSSI and MSSIM.   Figure  14) and GSSI (solid bar in Figure 14). The testing is done on a Dell computer with Intel(R) Core(TM) i7-4770 CPU, 16GB RAM (3.40 GHz). It is observed that the computational time of GSSI (3.92 seconds) is much lower than that of the MSSIM with lower order slider window size (2 Â 2 sliding window based MSSIM required 39.5 seconds). The results are as expected, the comparison of 20 Â 20-dimension OD matrices via 2 Â 2 sliding window needs to be performed (20 À 2 þ 1) Â (20 À 2 þ 1) ¼ 361 times. On the other hand, GSSI is an average value of all local SSI values computed 25 times (due to 25 geographical windows). In general, if we have G geographical windows and m x n sliding windows to compare M Â N OD matrices then, the GSSI approach is (Mm þ 1) Ã (N-n þ 1)/G times computationally more efficient than m Â n sliding window based MSSIM. From this, it can be concluded that GSSI is computationally more effective than MSSIM (based on sliding window of smaller sizes).

Sparse OD matrices and geographical windows from different zonal levels
The previous section demonstrated GSSI's application on dense OD matrices (SA2) and is based on geographical window from a particular higher zonal level (SA4). However, in most transport modeling applications we encounter OD matrices that are mostly sparse with zones much smaller than SA3. In such situations, it would be interesting to see how the GSSI approach performs for different types of zonal windows, and if it can compare sparse matrices better than the traditional MSSIM. Since BOD matrices constitute only partial observed flows, sufficient sparsity can be observed at SA2 level. Thus, the sparse OD matrices we considered for this analysis are B-OD matrices at SA-2 level (99 Â 99). A typical Monday (sparsity around 20%), is compared with typical Saturday and Sunday (sparsity around 16% and 21%) matrices. The different geographical window levels used are SA3 and SA4 as shown in Figure 15. The MSSIM considers only single window because the order of re-arranged OD matrices is found to be different during the three days. Refer to Figure 16 for the visualization of sparsity in SA3 (30101 to 30401, 30402, 30403, and 30404) and SA4 (301 to 304) windows. Table 3 demonstrates the GSSI comparisons applied using SA3 and SA4 windows as against the MSSIM values on SA2 OD matrices. We can see that MSSIM resulted in the highest similarity value but GSSI using SA4 and SA3 windows resulted in lower and lowest values of similarity, respectively. For example, Monday-Saturday comparison values are 0.9190, 0.8185, and 0.7202 by MSSIM, GSSI_SA4 and GSSI_SA3, respectively. The STR values (3 rd column in Table 3b) for Monday-Saturday (0.9730) and Monday-Sunday (0.9688) are so close to each other that there are hardly any structural differences between Sunday and Saturday travel patterns when compared to Monday. This is mainly due to the reason that the sensitivity of the SSIM formulation toward structural distortions reduces as the window size increases. Thus, we can see that GSSI is able to recognize the structural differences better than MSSIM despite sparsity in OD matrices.  Figure 13. Comparison of local travel patterns during Ekka, Saturday and Sunday.
Regarding the selection of geographical windows (i.e., consideration of either SA4 or SA3 as shown in Table 3), the smaller the size of window, the better is the ability of SSIM to compute local structural differences. Thus, we recommend using geographical windows from immediately higher level in the zonal hierarchy; for instance, use SA4 windows for SA3 OD matrices, SA3 for SA2 or SA2 for SA1, etc.

Sensitivity analysis
Although the mathematical formulation of GSSI is same as that of MSSIM, they vary in their similarity values due to different choice of local windows. Also, we have demonstrated in the previous section that GSSI is able to capture better structural distortions (as compared to MSSIM) among daily OD matrices. This section further tests the robustness of both GSSI and GSTR through sensitivity analysis for different conditions of OD matrix combinations. A statistical measure is considered to be robust if it is sensitive to changes in both "OD flows" and the "OD structure." The study site, data and the design of experimental set up are briefly discussed below.
1. Study site and data: The BCC region is the study area as shown in Figure 7. The reference OD matrix, X is the Bluetooth based OD matrix (20 Â 20) observed on Monday, March 7, 2016 (refer to Figure 8(a)). The query OD matrices are developed specific to the experiments and are obtained by perturbing the reference OD matrix. 2. Experiments: The experimental set up for the sensitivity analysis of GSSI and GSTR is shown in Table 4.
We design two sets of experiments for the sensitivity analysis. Each experiment tests the sensitivity of both GSSI (X, Y) and GSTR (X, Y) to capture the structural similarity and the structure of OD matrices, respectively. The experiments are as follows: a. Experiment based on uniform scaling effect: Here, the query matrices have the same structure as that of reference OD matrix while the magnitude of OD flows vary. If the uniform scale factor is one, then both OD matrices are exactly similar. Thus, abovementioned situation-1 and situation-2 are tested here (refer Table 4). b. Experiment based on random scaling effect: Here, the structure and OD flows vary between query and reference OD matrices. Thus, above-mentioned situation-3 is tested here (refer Table 4).
In the following sections, we first discuss the design of experimental set up for both experiments, and then present the results of sensitivity analysis.

Experimental set up for the sensitivity analysis of GSSI and GSTR
Criteria for uniform scaling effects Here, sensitivity of GSSI and GSTR are tested for different uniform scaling percentages. The reference OD matrix, X is compared with Y i where Y i ¼ u Ã X, and u is chosen from [0.1, 0.2, 0.3 … 1.9, 2.0]. GSSI is robust metric if: 1. The GSTR is zero for any value of u: 2. The GSSI value is zero for u ¼1 and should decrease as u deviated from unity.

Criteria for random scaling effects
Here, sensitivity of GSSI and GSTR are tested for four different cases of random scaling percentages i.e., w ¼ [5%, 10%, 15%, 20%] over three types of demand scenarios. These demand scenarios are generally encountered in traffic demand modeling (refer Djukic et al. (2015)) and are as follows: 1. Outdated surveys (low demand), 2. The best historical estimates (medium demand), and 3. Congested traffic conditions (high demand).
Note that the scenarios are named as low (l), medium (m), and high (h) in reference to the total daily demand on the network, and do not refer to the demands of individual OD pairs. In each case of the demand scenario, reference OD (X) is compared with 100 replications of query ODs (Y). The details of the demand scenarios are as follows: Low demand scenario: Here, both GSSI and GSTR compare X and where Y l i, w ¼ XÃð0:60 þ w Ã rand½0, 1Þ and i 2 ½1, 100: For instance, if w ¼20%, then Y l i, w ranges between 60% and 80% of X, and similarly for other values of w:

Results of uniform scaling effects
The results of uniform scaling for both GSSI and GSTR are shown in Figure 17. The plot illustrates that GSSI satisfies the conditions specified in the previous section i.e., GSSI values increased from 0.04 to 1 for 0.1 < ¼ u < 1 and decreased from 1 to 0.64 for 1 < u < ¼2. Similarly, GSTR remained unaffected (i.e., equal to 1) for both scaling-up and scaling-down cases (as described in the previous section). Thus, it is proved that both GSSI and GSTR are robust toward uniform scaling effects.

Results of random scaling effects
The plot shown in Figure 18(a) demonstrates decrease in the structural similarity as the magnitude of random fluctuations increase. For instance, GSSI for low demand scenario are 0.7759, 0.7675, 0.7489, and 0.7307 for w¼ 5%, 10%, 15%, 20%, respectively. A similar decreasing trend is observed for other demand scenarios. The structure component of GSSI, that is, GSTR also illustrates similar results in Figure 18(b). Both GSSI and GSTR are further compared using F-statistic to check if there is any difference among the means. The results of ANOVA shown in Table 5 demonstrate that F values are far greater than 3.24 at 0.05 level of significance. The p value for each F value is close to zero as shown within the brackets in Table  5. Therefore, the null hypothesis is rejected implying that the means of all four cases in each demand scenario are statistically different among each other for both GSSI and GSTR. Thus, the results prove that both GSSI and GSTR are robust toward random scaling effects.
Therefore, sensitivity analysis conducted using different OD matrix combinations in different situations should give enough confidence on the robustness of the proposed metric.

Discussion
The results of sensitivity analysis proved that the proposed GSSI approach is sensitive to both uniform and random demand fluctuations. Sensitivity toward uniform scaling effects prove that GSSI can be used to compare OD matrices of similar travel patterns (structures) such as comparing OD matrices of typical Monday and typical Tuesday. The experiments related to random scaling effects suggest that GSSI is able to capture random fluctuations at different demand scales (low, medium and high). The query matrices in high demand scenario are designed such a way that they are closer to reference matrix as compared to other scenarios. The values of GSSI can reflect proximity with higher GSSI (and GSTR) values in high demand scenario as compared to any other scenario. Sensitivity of GSSI toward random demand fluctuations demonstrate its ability in various applications such as (a) comparing target and estimated demands in OD estimation problem, (b) benchmarking OD estimation methods, (c) checking the convergence/ divergence of OD optimization algorithms, and (d) comparing the travel patterns based on type of the day (Behara et al., 2018) or modes (Hussain et al., 2019) . The ANOVA test results further proved that the mean values of GSSI and GSTR are statistically different for all four cases in each demand scenario. The study uses SA zones to demonstrate the concept. The SA zones are mainly defined considering the land use and population density. Travel behavior and demand is influenced by land use patterns and population density. Referring to the ATAP (2016) guidelines, "Transport zones should ideally contain homogeneous land use (for example, solely residential, industrial or commercial use or parking lots) and they should not cross significant barriers to travel (such as rivers, freeways and rail lines), and should have reasonably homogeneous access to the modeled transport systems. In this context, transport zones should match, as far as practically possible, Australian Bureau of Statistics (ABS) Statistical Area boundaries as defined by the Australian Statistical Geography Standard (ASGS). In general practice, transport zones  will be aggregations of the Statistical Area Level 1 (SA1) or mesh block boundaries." The Transport analysis zones (TAZ) used by Queensland Department of Transport and Main Roads (TMR) are well aligned with majority of SA zones and the household travel survey is designed at SA1 level; thus SA1 zones are linked to the Brisbane TAZs in the transport models (Weston, 2019). Although we considered SA zones as a proxy for TAZ, the concept of geographical windows still holds good should we have different levels of TAZs; for instance, consider TAZs used by the Florida Department of Transportation Systems Planning Office (2007). It uses three types of TAZs for three different purposes: large-sized zone for state-wide planning, medium-sized zone for arterial planning, and small-sized zone for corridor analysis. To apply the GSSI's rationale, medium sized TAZs (with approx. 4000 persons per zone) can be used to define geographical windows for small-sized TAZs (with approx. 1200 to 3000 persons per zone).

Conclusion
A holistic comparison of OD matrices should include both deviations of individual OD flows and the structural (dis)similarity between the OD matrices. The traditional statistical measures compare OD matrices using deviations of individual OD flows, and fail to account for OD matrix structural differences. Very limited studies in the past focused on developing such structural (dis)similarity measures (refer to Table A1 in the Appendix section for conceptual comparison of these measures). These metrics are popular in other disciplines and are introduced in transport applications. One such metric is MSSIM that computes statistics on groups of OD pairs defined by local windows. However, comparison of OD matrices using MSSIM can result in different values based on the choice of local window size. In literature, no clear consensus has been reported on the level of acceptability of sliding window size and the resulting MSSIM values. While a few researchers proposed variants of MSSIM, they needed further investigation. We propose an alternative way of addressing the limitations of MSSIM with respect to local window size and geographical adjacency through the concept of geographical window defined by geographical boundaries. We term this extension to MSSIM as geographical window based structural similarity index (GSSI). The findings of this study are (a) GSSI technique can capture local travel patterns better than a sliding window based MSSIM, (b) it has a physical meaning to local statistics, (c) it is computationally more effective, (d) it could compare large-scale sparse OD matrices, and (e) it can identify subtle structural differences among OD matrices. The study suggests that it is better to choose geographical windows at a level immediately higher to the zonal level of OD matrices that are to be compared. A thorough sensitivity analysis proves that GSSI is a robust statistical measure.
As a part of future research directions, we plan to study the optimality of stability constants in the SSIM formulation and weights specific to SSIM components (this is important while comparing sparse and dense OD matrices). Although, we have tried to address the issue related to selection of window size, the level of acceptability of GSSI still depends on different transport applications, and this requires further discussions with practitioners and scientific community in the transport domain.
The proposed formulation requires knowledge of the network in terms of its statistical boundaries and its hierarchy. Such information is generally prepared and maintained by the Bureau of Statistics and should be available to the transport analyst. In the current application, the hierarchical structure of the zones provided by the Australia Bureau of Statistics is utilized to define the dynamic window size. For future study, it is recommended to establish guidelines for the definition of the geographical boundaries. It should consider zonal socio-economic characteristics, land use pattern, and spatial proximity.