{'fraction_non_alphanumeric': 0.02921876079298197, 'fraction_numerical': 0.00511155626165642, 'mean_word_length': 4.9103898754847926, 'pattern_counts': {'":': 0, '<': 1, '<?xml version=': 0, '>': 2, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, '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}}
From Collective Adaptive Systems to Human Centric Computation and Back: Spatial Model Checking for Medical Imaging

Recent research on formal verification for Collective Adaptive Systems (CAS) pushed advancements in spatial and spatio-temporal model checking, and as a side result provided novel image analysis methodologies, rooted in logical methods for topological spaces. Medical Imaging (MI) is a field where such technologies show potential for ground-breaking innovation. In this position paper, we present a preliminary investigation centred on applications of spatial model checking to MI. The focus is shifted from pure logics to a mixture of logical, statistical and algorithmic approaches, driven by the logical nature intrinsic to the specification of the properties of interest in the field. As a result, novel operators are introduced, that could as well be brought back to the setting of CAS.


Introduction
Formal verification of properties of Collective Adaptive Systems (CAS) is a challenging subject. The huge number of considered entities introduces a gap with classical finite-state methods as the number of states grows exponentially. Approximation methods, such as mean-field or fluid-flow approximation, have been proposed to mitigate this aspect (see [5,4,30]). Another relevant issue is that of spatial distribution of the considered entities. Entities composing a CAS are typically located, and moving, in a physical or logical space. Collective behaviour is driven by interaction, which is frequently based on proximity. This makes spatial aspects more prominent in the case of CAS than in classical concurrent systems, and leads to the introduction of spatial properties in formal verification.
Model checking [1] is a formal verification technique that is based on static analysis of system properties that are described by modal logics. Modal operators have been traditionally used to denote possibility or necessity, probability, constraints on continuous time, access to security contexts, separation of parallel components of a system and so on. But since the very beginning, modal logics have also been interpreted on spatial structures, such as topological spaces (see [3] for a thorough introduction). In this context, formulas are interpreted as sets of points of a topological space, and in particular φ is usually interpreted as the points that lay in the closure of the interpretation of φ . A standard reference is the Handbook of Spatial Logics [33]. Therein, several logics are described, with applications far beyond topological spaces; such logics treat not only aspects of morphology, geometry, distance, but also advanced topics such as dynamic systems, and discrete structures, that are particularly difficult to deal with from a topological perspective.
Model checking of spatial (and spatio-temporal) logics is a more recent development (see e.g., [15], [22], [23]). In [12], a model checking algorithm for spatial logics in the topological sense has been proposed. Therein, models are based on extensions of topological spaces called closure spaces, designed to accommodate also discrete structures such as finite graphs in the topological framework. The so-called surrounded operator is introduced to describe points that are located in a region of points satisfying a certain formula and surrounded by points satisfying another one (in other words, the surrounded operator is a spatial form of the until operator of temporal logics). The resulting logic is called SLCS (Spatial Logic of Closure Spaces). In [11], a spatio-temporal model checking algorithm using the same principles has been proposed, for a logic combining the spatial operators of SLCS with the well-known temporal operators of the Computation Tree Logic CTL. Applications to CAS have been developed in the fields of smart transportation [10], and bike-sharing systems [13]. A free and open source implementation of the spatio-temporal model checking algorithm of [11] is provided by the tool topochecker 1 .
In this position paper, we explore various ideas for the application of spatial and spatio-temporal model checking that depart from the setting of CAS and are tailored to Human Centric Computation, in particular to the field of Medical Imaging (MI). In this domain, space consists of points called voxels, that are arranged in a multi-dimensional, possibly anisotropic grid. Several spatial analyses are performed in MI. Such analyses are typically described by structured combinations of attributes related to proximity, shape, aspect and distance of features of interest, for which spatial logics provide a well-suited descriptive language. Our work is motivated by some considerations about how medical image analysis is carried out. The overall general description of a feature (e.g. the shape and spatial arrangement of parts of an image that exhibit diseases) is often carried out informally, but in a logically structured way (e.g.: "the tumour is lighter than the surrounding brain area, and touches the oedema, whose intensity is a bit darker than the tumour"). Such description is then turned into a series of different analysis passes, sometimes performed by specific tools, but frequently done by writing ad-hoc programs. The results of such different passes are often integrated by hand or using hand-crafted scripts. This complex and elaborate process hinders the implementation, and sharing across the medical community, of novel analysis methods that emerge from current research. A leap forward is needed in the unambiguous and precise specification of such procedures, which could be provided by logical methods borrowed from Computer Science and in particular the area of formal methods. Our research program aims at paving the way and establishing foundational results for such a development to happen.

Spatial logics for medical imaging
The contribution of Computer Science to the field of medical image analysis is increasingly significant, and will play a key role in future healthcare. Computational methods are currently in use for several different purposes, such as: Computer-Aided Diagnosis (CAD), aiming at the classification of areas in images, based on the presence of signs of specific diseases [17]; Image Segmentation, tailored to identify areas that exhibit specific features or functions (e.g. organs or sub-structures) [20]; Automatic contouring of Organs at Risk (OAR) or target volumes (TV) for radiotherapy applications [6]; Indicators finding, that is, the identification of indicators, computed from the acquired images, that permit early diagnosis, or understanding of microscopic characteristics of specific diseases, or help in the identification of prognostic factors to predict a treatment output [9,43]. Examples of indicators are the Mean Diffusivity and the Fractional Anisotropy obtained from Magnetic Resonance (MR) Diffusion-Weighted Images, or Magnetization Transfer Ratio maps obtained from a Magnetization Transfer acquisition [16,31].
Such kinds of analyses are strictly related to spatial and temporal features of acquired images. For example, the diagnosis of specific diseases requires the observation of variations in time of the response to a particular acquisition technique. Another example is provided by longitudinal studies, that consist of repeated observations of the same variables over long periods of time, to help understanding the areas involved in some diseases. Such investigations can take advantage of spatial and temporal capabilities in model checking. In this work, we provide a preliminary investigation of spatial logical operators that may be used to identify areas in the images based on local or non-local features and prior knowledge (e.g., areas surrounded by or near to or similar to areas with particular properties).
We start from the spatial logic SLCS presented in [12]. SLCS features the spatial operators near and surrounded. The syntax of SLCS is defined by the following grammar, where p ranges over P, namely the set of atomic propositions: Φ :: In the syntax, is the constant true; ¬ and ∧ are standard logical negation and conjunction. Formula N φ denotes all points that are "near" to the interpretation of φ , whereas formula φ 1 S φ 2 denotes all points that lie in a subset of the interpretation of φ 1 , whose "boundary" satisfies φ 2 . The notions of "near" and "boundary" are made formal in the interpretation of the logic, resorting to a closure model [12] for details), and v : P → P(X) is the valuation of atomic propositions.
In the remainder of this paper, we discuss some additional logical operators that deal with distance (Section 3) and texture analysis (Section 4), and we reconsider the models of SLCS in the light of such additions. These operators are implemented in the experimental branch of topochecker. In models associated to medical images, X is the set of points of an image, and C associates to each subset A of X the union of the neighbours of each point of A, by a user-defined notion of neighbourhood. To cope with the quantitative information present in medical images, atomic propositions have quantitative valuations over the real numbers, instead of just boolean valuations. However, to retain the boolean interpretation of SLCS formulas, the syntax of atomic propositions in the logic is that of constraints with variables ranging over R. Atomic predicate p is a shorthand for p = 1. For example, formula p > k ∧ q < h denotes all points x ∈ X such that the value of p is greater than k and the value of q is less than h. Formally, this does not require changes to SLCS and its semantics: given a setP of proposition letters, with quantitative valuationv :P → R X , the set of atomic propositions P is given by the set of all possible constraints onP, and the valuation function v : P → 2 X is just evaluation of constraints, which makes use ofv.

Distance operators
Distance operators can be added to spatial logics in various ways (see [28] for an introduction). Distances are very often expressed using the real numbers R. Typically, one considers operators of the form D e(z) φ where e(z) is a constraint parametrised by a free variable z, and φ is a formula denoting a spatial property. The intended semantics is that point x is a model of D e(z) φ if and only if there is a point y satisfying φ such that the distance d from x to y satisfies the constraint e(d). Logics of metric spaces have been introduced in [29]; therein, the constraint e(z) can only be in the form z ≤ k, where k ∈ R (distance at most k), or k 1 ≤ z ≤ k 2 (distance included between k 1 and k 2 ) 2 . The latter is called "doughnut operator" in [28]. Notably, the doughnut operator cannot be expressed just using z ≤ k in combination with boolean operators. In the following, we will discuss a general model checking procedure which is able to verify the satisfaction of such formulas in an efficient way.
Models. First, we need to discuss appropriate models for our logics. The quasi-discrete closure models of [12] provide a starting point. A quasi-discrete closure model is completely described by a pair (X, R) where R is a binary relation on a set of points X -that is, a (possibly directed, unlabelled) graph -and by a valuation v associating to each atomic proposition, in a finite set P, a set of points of X. This data uniquely defines a model M = ((X, C ), v) where C : 2 X → 2 X is the closure operator, that coincides with the dilation operation of the graph (X, R). Such models do not include information on distances. One possibility is to enrich the structure of M by adding a distance operator d : X × X → R ≥0 . In medical imaging, the structure of a metric space is a very natural setting, as typical distances are based on either Euclidean spaces or symmetric graphs. See Appendix A for more details on the possible choices for a specific metric in medical imaging applications. An open question is what are the additional axioms (one might say "compatibility conditions") linking closure and distance. The link is well-known for the case of topological spaces. More precisely, topological spaces are obtained from metric spaces by defining open sets as those sets S such that all points of S have a neighbourhood 3 in S; a generalization to closure spaces is an interesting topic for future research.
Distances. In the case of quasi-discrete closure spaces generated by a graph (X, R), it is natural to consider distance operators that are obtained by the shortest path distance of a weighted graph obtained by assigning weights to the arcs in R. In [36], such models are used for spatial logics featuring a metric variant of the surrounded operator of [12]. Note however that in some cases other notions of distance can be more appropriate. For example, sampling an Euclidean space is often done using a regular grid, where points of a graph are arranged on multiples of a chosen unit interval, and connected by edges using some notion of connectivity (e.g. in 2-dimensional space, one typically uses four or eight neighbours per point). Shortest-path distance and Euclidean distance obviously divert in this case (see Figure 2 and Figure 3 in Appendix A), no matter how fine is the grid or how many finite neighbours are chosen. Medical images introduce some more complexity due to the fact that multi-dimensional voxels frequently are anisotropic, that is, the sizes of a voxel in each of the multiple dimensions of the image are different.
Global model checking through distance transforms. The model checking algorithm implemented in topochecker is a global one. That is, all points of the considered space or space-time are examined, and those satisfying the considered formula are marked by the algorithm. Even though this classical approach has drawbacks -mostly related to the restriction to finite models, and the fact that large parts of a model need to be stored in central memory -it is very helpful in the case of medical image analysis, as logical operators may take advantage of global analysis of the whole space. One example where this is particularly useful is the computation of distance formulas. This can be done using so-called distance transforms. The concept of distance transform comes from the area of topology and geometry in computer vision [27]. The idea is extensively used in modern image processing. Given a multi-dimensional image equipped with Euclidean distance, global spatial model checking of formulas D e(z) φ can be done in linear time with respect to the number of points of the space, assuming that the computation time of e(z) is negligible, and that the computation of φ is linear in turn, which is true for the spatial logic SLCS of [12]. Consider a multi-dimensional image. The outcome S of computing the truth value of φ on each point of a model is a binary (multi-dimensional) image. The binary value stored in each point corresponds to the truth value of φ on that point. From this data, it is possible to define a transformed image, called the distance transform of S, such that in every point x, a value d x ∈ R is stored. The value of d x is meant to correspond to the minimum distance between x and a point satisfying φ . There exist both exact and approximate algorithms for computing distance transforms, that are usually classified by their asymptotic complexity, computational efficiency, and possibility of parallel execution. In particular, there are effective linear-time algorithms [34,18]. In global model checking, one first computes the distance transform, and then in one pass, for each x, the quantitative value of d x can be replaced with the boolean value of e(d x ). It is worth noting that similar algorithms exist for weighted graphs equipped with shortest-path distance. In this case, asymptotic complexity is generally speaking quasi-linear but not linear, although execution time is highly dependent on the structure of the considered graph.
In the next section we shall discuss how distance operators can be combined with texture analysis to identify tissues of different nature that lay within a certain distance of each other.

Texture analysis operators
Texture analysis (TA) operators are designed for finding and analysing patterns in medical images, including some that are imperceptible to the human visual system. Patterns in images are entities characterised by brightness, colour, shape, size, etc. TA includes several techniques and has proved promising in a large number of applications in the field of medical imaging [26,32,7,14]; in particular it has been used in CAD applications [44,24,25] and for classification or segmentation of tissues or organs [8,39,38]. In TA, image textures are usually characterised by estimating some descriptors in terms of quantitative features. Typically, such features fall into three general categories: syntactic, statistical, and spectral [26]. Our preliminary experiments have been mostly focused on statistical approaches to texture analysis. Statistical methods consist of extracting a set of statistics descriptors from the distributions of local features at each voxel. In particular, we studied first order statistical methods, that are statistics based on the probability density function of the intensity values of the voxels of parts, or the whole, of an image, approximated as a histogram collecting such values into batches (driven by ranges). In this approach, the specific pixel adjacency relationship is not taken into account. Common features are statistical indicators such as mean, variance, skewness, kurtosis, entropy [40]. Although a limitation of these operators is that they ignore the relative spatial placement of voxels, statistical operators are important for MI as their application is invariant under transformations of the image. In particular, first order operators are, by construction, invariant under affine transformations (rotation and scaling), which is necessary when analysing several images acquired in different conditions. Nevertheless it is possible to construct features using first order operators, keeping some spatial coherence but losing at least partially the aforementioned invariance [42].
In our experimental evaluation, we defined a logical operator, called SCMP -for statistical comparison -that compares areas of an image that are statistically similar to a predetermined area SA, identified using a sub-formula. More precisely, SCMP is used to search for sub-areas in the image whose empirical distribution is similar to that of SA, up-to a user-specified threshold. For each voxel, a small surrounding area is considered; its statistical distribution is compared to that of SA and a threshold is applied, obtaining a Boolean value that denotes whether the voxel belongs to an area statistically similar to SA.
Statistical distributions are compared using the cross-correlation function. Currently, the syntax and interpretation of the operator are quite experimental, and different methods may be used to define the data aggregation and comparison of distributions which are needed for its implementation. Complexity is, among others, an issue to deal with.
The SCMP operator results in a generalisation of classical TA based on first-order statistics, since it analyses the statistical distribution of a neighbourhood of each voxel as a whole, whereas classical techniques for TA resort to the extraction of specific indicators. In Figure 1 we show the output of topochecker, enhanced with the SCMP operator, and the distance-based operators described in Section 3, applied to a slice of a MR-Flair acquisition of a brain affected by glioblastoma (GBM), an intracranial neoplasm 4 . GBMs are tumors composed of typically poorly-marginated, diffusely infiltrating necrotic masses. Even if the tumor is totally resected, it usually recurs, either near the original site, or at more distant locations within the brain. GBMs are localised to the cerebral hemispheres and grow quickly to various sizes, from only a few centimetres, to lesions that cover a whole hemisphere. Infiltration beyond the visible tumor margin is always present. In MR T2/Flair images GBMs appear hyperintense and surrounded by vasogenic oedema 5 .
Being able to segment tumor and oedema in medical images can be of immediate use for automatic contouring applications in radiotherapy and, in perspective, it can be helpful in detecting the invisible infiltrations in CAD applications. In Figure 1b we show a segmentation obtained using topochecker. First, two thresholds are applied to the original image, in order to identify two areas of the image with particular brightness, that are supposed to loosely correspond to an oedema and a tumor. Furthermore, the distance operator (shortest path distance with 9-neighbours, see Appendix A) is used to imposes a certain degree of proximity between the oedema and the tumor. The voxels assigned to the oedema are drawn in yellow; the voxels assigned to tumor are drawn in orange. These areas are used as two different SA for the analysis shown in Figure 1c, obtained using the SCMP operator to further enlarge the two regions, by searching areas that, although not falling in the specified thresholds, are still identifiable as oedema and tumor. The model checker colours the additional voxels assigned to tumor in blue; these areas are characterised by a statistical distribution similar to the orange region. The tool colours the additional voxels assigned to oedema in magenta. These are obtained by searching for areas having statistical distribution similar to the yellow part of the image. In Figure 1d we show the final result of the procedure with tumor voxels in orange and oedema voxels in yellow.

Discussion
We have just started the exploration of logical methods for medical image analysis in the domain of radiotherapy. Logical properties are used as classifiers for points of an image; this can be used both for colouring regions that may be similar to diseased tissues, and therefore being diseased tissue in turn, and for colouring regions corresponding to organs of the human body. Envisaged applications range from contouring to computer-aided diagnosis. The field of Spatial Logics can benefit from such kinds of research; for example, texture analysis operators may be defined, and in particular operators that compare regions based on statistical similarity. It will be interesting to study relevant axioms and theories that can deal with the uncertainty generated by analysis of statistical similarity in theorem proving or completeness studies, as well as theories of bisimilarity or minimisation of models. From the model-(a) A slice of a FLAIR MR acquisition of a brain affected by a glioblastoma.
(b) Same slice superimposed to regions of interest (ROIs) determined by threshold and distance operator. The yellow ROI is the oedema; the orange one is the tumor.
(c) Same slice superimposed to ROIs determined by the SCMP operator, starting from orange and yellow ROIs. In blue, the additional tumor voxels; in magenta, the additional oedema voxels.
(d) Final result of application of threshold, distance, and the SCMP operator.
The yellow area is the identified oedema; the orange area is the tumor. theoretical point of view, questions of interest relate to the various kinds of distances that may arise in the considered spatial models, ranging from classical Euclidean distance to shortest-path distance on weighted graphs, on their axioms, and the relation between different notions. The effect of approximate distance transforms on the representational complexity of models may be an interesting question for future work. Our early experiments show that typical analyses carried out using spatial model checking in medical imaging require careful calibration of numeric parameters (for example, a threshold for the distance between a tumor and the associated oedema, or the size of areas identified by a formula, that are small enough to be considered noise, and ought be filtered out). The calibration of such parameters can be done using machine-learning techniques. In this respect, future work could be focused on application, in the context of our research line, of the methodology developed by Bartocci et al. (see e.g. [19,2]). Some recent research focused on themes that are close to our planned development. In particular, [41] uses spatio-temporal model checking techniques inspired by [22] -pursuing machine learning of the logical structure of image features -to the detection of tumors. In contrast, our approach is more focused on human-intelligible logical descriptions. The work [37] is closer to the setting of CAS, applied to biological processes, with an interesting focus on multi-scale aspects. The research line that we present in this paper stems from research in collective adaptive systems and departs from it to direct spatial analysis to medical imaging. We foresee that the novel statistical texture analysis operators and the study of global model checking of distance formulas using distance transforms are of interest when dealing with very large populations that are spread over some spatial structure (e.g. a geographical map). Potential applications include the analysis of statistical properties arising from gossip protocols and disease spreading models, in which statistical distribution of features in space appears to be relevant.

A Comparison of metrics
In MI, Euclidean distance (based on the 2-norm) is the reference distance between two voxels. Therefore, in this context, Euclidean distance is considered error-free. We shall now discuss different alternative definitions of distance that could be interesting for spatial model checking, both for medical imaging applications and also as a valuable addition to spatial analysis methods for CAS. Starting from [21], we experimented with two distance transform operators, the multi-dimensional Euclidean error-free distance operator, called EUCL in the following, and a distance transform operator based on Dijkstra's shortest path algorithm, called MDDT (Modified Dijkstra Distance Transform), that operates on the shortest-path distance of a graph constructed from the image. From the point of view of the Dijkstra algorithm, an image is a graph whose vertices are the voxels and whose arcs connect each voxel x with chosen voxels that are considered adjacent to x. Every arc is labelled with the chosen distance function applied to the two vertices that the arc connects. In other words, we are considering a graph with nodes located in a distance space, arcs weighted according to the distance of the space, and a chosen notion of adjacency. In this particular case, the shortest-path distance is also called Chamfer distance. The chosen adjacency is the most important factor in the precision-efficiency trade-off of the computed distance: the more adjacent voxels are considered, the more precise is the Chamfer distance when compared to the Euclidean distance, at the expenses of generating graphs with larger out-degrees. In Figure 2 and Figure 3 we show in red the points satisfying D z>k φ for a binary image with only one point satisfying φ (in the centre of the image). In Figure 2 we show the output of the error-free EUCL operator. In Figure 3, we show the output of the MDDT operator, alongside the characteristic pattern of the percentage error with respect to Euclidean distance (Figure 3b and 3d -see [21] for a detailed analysis of the percentage error of several distance transform algorithms). The percentage error δ (x) for the distance transform d(x) is defined in every voxel x as δ (x) = |d eucl (x)−d(x)| d eucl (x) . In Figure 3a and Figure 3c we show the Chamfer distance obtained using MDDT, with different choices of adjacent voxels. In Figure 3a the adjacent voxels of a point are chosen to be its immediate neighbours on the main directions and diagonals (called Moore neighbourhood in 2d images). In Figure 3c, again the main directions and diagonals are used, but this process is iterated two times (that is, the chosen adjacent voxels of x are all points in a hypercube of size 5, centred on x, except x itself).
In a Euclidean space, Euclidean distance is not the only possible distance; several definitions exist based on different norms. In (Figure 4) we depict two widely used metrics. In Figure 4a adjacency is the same as in Figure 3a but all arcs have weight 1. This is called the Chebyshev or chessboard distance, that in a Euclidean space is the distance based on the so-called infinity-norm. In Figure 4b, the underlying graph is as in Figure 4a but adjacency contains only voxels whose coordinates differ at most in one dimension) (for 2D images, this is called Von Neumann neighbourhood), thus obtaining the taxicab or cityblock distance. In a Euclidean space, the cityblock distance is the distance based on 1-norm.