"Saliency maps for of images" Cross-Reference to Related Applications [0001] The present application claims priority from Australian Provisional Patent Application No 2022902976 filed on 11 October 2022, the contents of which are incorporated herein by reference in their entirety. Technical Field [0002] This disclosure relates to creating saliency maps for the classification of images. Background [0003] Machine learning models, such as deep neural networks, are powerful tools in addressing challenging medical imaging problems. For example, machine learning models can be used to determine the presence of a disease from an image or a collection of images. However, for an ever-growing range of applications, interpreting a model’s prediction remains non-trivial. For example, while the models may be able to detect the presence of a disease from image data, it may not be clear how the model makes its prediction nor may it be clear what features the model deems important when making its prediction. The understanding of decisions made by black-box algorithms is a challenge and assessing their fairness and unbiasedness is an important step to healthcare deployment. Therefore, what is needed is a method that can accurately determine how machine learning models make their prediction from image data. [0004] Any discussion of documents, acts, materials, devices, articles or the like which has been included in the present specification is not to be taken as an admission that any or all of these matters form part of the prior art base or were common general knowledge in the field relevant to the present disclosure as it existed before the priority date of each of the appended claims. [0005] Throughout this specification word “comprise”, or variations such as “comprises” or “comprising”, will be understood to imply the inclusion of a stated element, integer or step, or group of elements, integers or steps, but not the exclusion of any other element, integer or step, or group of elements, integers or steps. Summary [0006] A method for creating a saliency map associated with a first classification of image data by a classifier, to facilitate visual interpretation of the first classification by a user, the method comprising: applying the classifier to the image data to determine the first classification; distorting the image data to produce distorted image data, wherein applying the classifier to the distorted image data produces a second classification that is opposite to the first classification of the image data; and producing the saliency map based on the image data and the distorted image data, the saliency map indicating salient features of the image data in relation to the classifier. [0007] It is an advantage to distort the image data to produce a classification that is opposite to the first classification, as this gives the greatest contrast of the classifier’s classification. Therefore, producing a saliency map based on the image data and the distorted image data shows salient features that the classifier considers important when making the classification of the image data. [0008] In some embodiments, producing the saliency map based on the image data and the distorted image data comprises calculating a difference between the image data and the distorted image data. [0009] In some embodiments, the method further comprises determining a distortion of the image data by adjusting the distortion to optimise an output of the classifier applied to the distorted image data. [0010] In some embodiments, the further comprises repeating the optimisation until the output of the classifier indicates that the second classification is opposite to the first classification. [0011] In some embodiments, the distortion is parameterised by a set of grid points that span the image data and determining the distortion comprises iteratively updating the set of grid points until the output indicates that the second classification is opposite to the first classification. [0012] In some embodiments, optimising the output of the classifier further comprises optimising a cost function and calculating a gradient to update the grid points in the direction of the gradient, wherein the cost function is a function of the distortion. [0013] In some embodiments, distorting the image data comprises a free-form deformation. [0014] In some embodiments, distorting the image data with the free-form deformation comprises applying a deformation vector field to voxels of the image data. [0015] In some embodiments, determining the distortion comprises at least one penalty term to enforce a constraint on the distortion. [0016] In some embodiments, at least one penalty term is indicative of any one of: a local invertibility of the distortion; a sparsity of the distortion; and a subset of the image data. [0017] In some embodiments, the distorted image data is indicative of image data on a decision boundary of the classifier. [0018] In some embodiments, the classifier has exactly two possible output classifications. [0019] In some embodiments, the data comprises images of an anatomy of a patient. [0020] In some embodiments, the distorted image data is anatomically plausible. [0021] In some embodiments, the images form a spatial sequence of two-dimensional cross sectional images. [0022] In some embodiments, the spatial sequence of two-dimensional cross sectional images forms an input to the classifier to determine the first or second classification. [0023] In some embodiments, the image data comprises magnetic resonance image data or computed tomography image data. [0024] In some embodiments, the classifier is a neural network. [0025] Software that, when installed on a computer and executed by the computer, causes the computer to perform the above method. [0026] A system for creating a saliency map associated with a first classification of image data by a classifier, to facilitate visual interpretation of the first classification by a user, the system comprising a processor configured to: apply the classifier to the image data to determine the first classification; distort the image data to produce distorted image data, wherein applying the classifier to the distorted image data produces a second classification that is opposite to the first classification of the image data; and produce the saliency map based on the image data and the distorted image data, the saliency map indicating salient features of the image data in relation to the classifier. [0027] Optional features provided in relation to the method, equally apply as optional features to the software and the system. Brief Description of Drawings [0028] An example will now be described with reference to the following drawings: [0029] Fig.1 illustrates a method for creating a saliency map. [0030] Fig.2 illustrate an example of a spatial-sequence of two-dimensional cross-sectional images indicative of a MRI scan or CT scan. [0031] Fig.3 illustrates a system for creating a saliency map. [0032] Fig.4 illustrates a convolution applied to the data structure from Fig.2. [0033] Fig.5 illustrates a fully connected output layer for the classifier with an output value corresponding to the first or second classification. [0034] Fig.6a illustrates an overview of the disclosed method. [0035] Fig.6b illustrates the original image x , counterfactual image x ^ and the corresponding saliency map x ^ ^ x . [0036] Fig.7 illustrates the architecture of the DNN used to evaluate the disclosed method. [0037] Fig.8 illustrates robustness assessment of the disclosed method for Gaussian noise with a standard deviation of 0.1 injected on the input image. [0038] Fig.9 illustrates the saliency map x ^ ^ x overlayed on the on the initial patient’s MR Image. This image depicts the slight change of the initial image that results in the model changing its prediction. [0039] Fig.10 illustrates saliency maps for each method of an axial and a sagittal view given a single prediction. [0040] Fig.11 illustrates the success in changing the network’s prediction when extremal deformations are applied on individual substructures of the brain. Description of Embodiments [0041] During this last decade, medical Deep Learning (DL) models have surpassed deterministic approaches by delivering faster runtimes and above human-expert performance. Despite this, DL-based methods suffer slow adoption rates in healthcare settings, with critics often referring to their lack of tractability or explainability. To this end, numerous methods have been devised to provide the user with an understanding of the “why" by highlighting the most critical parts of the input contributing to the prediction. However, these methods often produce explanations that are difficult to interpret by a human user or lack robustness. [0042] The tractability or explainability of DL-based models could be also used for medical research purposes. One example may be that a DL-based models may be able to detect a disease from a patient image that would otherwise be undiagnosed by a medical professional. Understanding how the DL-based methods made its diagnosis prediction may lead to a better understand of the diseases, such as the discovery of new symptoms associated with the disease. This could lead to better diagnosis and treatment of the disease in the future. [0043] One way to determine the tractability or explainability of machine learning models, such as DL-based methods, is through the use of saliency maps. Saliency refers to the importance of a particular feature in an image for an output classification, as determined by the machine learning models. In that sense, a saliency map is a two- dimensional representation of the influence that each point or area in the map has on the output classification. More specifically, in the area of image processing, a saliency map is a two-dimensional image where each pixel has a pixel value (which may be represented as a brightness or colour) according to the importance of that pixel to the output classification. [0044] For example, the saliency map be represented by a greyscale image where each pixel has an intensity value. In this example, a pixel with a high intensity value has a large influence on the output classification while a pixel with a low intensity value has a low influence on the output classification. The saliency map may also be applied as an image mask to the original image to display only areas of the original image that influence the classification. Equivalently, a saliency map can be a three- dimensional representation of the influence that each point or area in the map has on the output classification, in such an example where the machine learning models produce a classification of a three-dimensional image. [0045] Some existing methods are computationally intensive, especially for three- dimensional images such as computed tomography (CT) images or magnetic resonance (MR) images, and can be challenging to train, particularly on small datasets. [0046] In this disclosure, a deformation based counterfactual explanability method is presented that does not possess the inherent limitations of other methods. The disclosed method comprises an extremal perturbation method which, given an input image and a classifier (machine learning model), such as a black-box Deep Neural Network (DNN) for example, produces a delusive image example by optimising a regular deformation or distortion on the input image. The image is distorted in such a way that when the classifier is applied to the distorted image, it produces the opposite classification obtained from the original image. [0047] From the deformation or distortion of the image, a saliency map (pixel attribution map) is derived, which indicates the salient (important) features that the machine learning model considers in the image when classifying the image and making a prediction about the image. The contrast between the original image that has an associated first classification and the distorted image that has an associated second classification, where the second classification is opposite to the first classification, allows salient features of the original image to be accurately determined. This is because the resultant saliency map visually shows the pixels of the original image that are used in differentiating the first classification from the second classification. [0048] These pixels are often spatial or clustered together to form a salient area or region, where this area or region of pixels may be referred to as a salient feature. It is possible that there are one or more salient feature in a single saliency map, which each may influence the classification determined by the classifier. These salient features are then presented visually to a user, as a way to visually present the pixels which influence the output of the classifier, and are used to visually interpret the important areas of the image. [0049] Here, images may be represented as image data. Image data may be indicative of a single image or multiple images. The multiple images may also have an associated order or sequence. [0050] The disclosed method is compared against classical post-hoc explainability methods, in order to determine the ability of the disclosed method to identify salient features. This disclosure also presents a more instructive visualisation technique based on the difference between the original input and the generated distorted image. The disclosed method is used to demonstrate the robustness of lightweight deep Convolutional Neural Networks (CNNs) against localised alterations of the input image. These results highlight the robustness of this architecture. [0051] The disclosed method is presented and tested through the use of a DNN trained to classify neurodegenerative diseases from brain MRI. It is noted that while the disclosed method is presented and tested using this DNN, the disclosed method is not limited to DNNs only. In this disclosure, a diffeomorphic transformation is optimised that deforms a given input image to change the prediction of the model. This provides anatomically meaningful, and thus easily interpretable by clinicians, saliency maps showing tissue atrophy and expansion that would change the model decision, in this case from healthy to Alzheimer’s disease (AD). This disclosed method is then compared to other saliency methods and shows that it provides more meaningful saliency maps when applied to one subject, but also disease consistent atrophy patterns when used over larger cohorts. [0052] The disclosed method does not domain knowledge. Here, domain knowledge refers to the existing body of knowledge (mathematical formulations) about the MR image reconstruction based on physical models such as parallel imaging. The disclosed method is also model agnostic, meaning that machine learning model is not limited to any particular model, such as the deep neural network trained to classify neurodegenerative diseases from brain MRI. [0053] A counterfactual example (distorted version of the original image) is produced through the resolution of an optimisation problem over a set of small and regular deformations of the input image. This disclosure evaluates this idea using 3D volumetric structural MR images, with the following assumptions: ^ The disease continuum can be modelled by smooth and invertible mappings. ^ The mapping from the healthy control to the disease group can be generated using a sufficiently refined Free-Form Deformation. [0054] Similar assumptions may be applied if the classifier related to applications outside of the medical field. For example, if the classifier was a multi-layered perceptron that classifies images of cats and dogs, it may be assumed that a Free-Form Deformation can be applied to produce the corresponding saliency map that can be modelled by smooth and invertible mappings. [0055] In this sense, distorting the image data comprises a Free-Form Deformation (FFD). FFD is a geometric technique used to model simple deformations and is based on the idea of enclosing an object within a cube or another hull object, and transforming the object within the hull as the hull is deformed. Deformation of the hull is based on the concept of so-called hyper-patches, which are three-dimensional analogs of parametric curves such as Bézier curves, B-splines, or NURBs. In essence, the pixel or voxels values of the image are perturbed slightly according to a fixed set of rules to achieve the desired deformation. [0056] Fig.1 illustrates a method 100 creating a saliency map. The saliency map is associated with a first classification of image data by a classifier, such as a trained machine learning model. The saliency map is used to facilitate visual interpretation of the first classification by a user. In method 100, a processor of a computer system first applies 101 the classification to the image data to determine the first classification. For example, the image data may correspond to a MR image of a patient’s brain and the classifier may classify the MR image as having Alzheimer’s disease or not having Alzheimer’s disease. In this example, the first classification may correspond to having Alzheimer’s disease. [0057] The processor then distorts 102 the image data to produce distorted image data, in such a way that applying the classifier to the distorted image data produces a second classification that is opposite to the first classification of the image data. In the example above, if the first classification corresponds to having Alzheimer’s disease, then the second classification would correspond to not having Alzheimer’s disease, thereby being opposite to the first classification. The processor may distort 102 the image data only slightly, so that there is a minimal change between the original image data and the distorted image data. The change may be minimal, such that, the original image and the distorted image are indistinguishable from one another by visual inspection. This gives an indication of the minimum requirements needed to change the classification of the classifier. In other words, the minimal change to the image data required to change the first classification to the second classification. Distorting 102 the image data to produce distorted image data may involve applying a distortion to the image data to produce the distorted image. The distortion may comprise a FFD, for example. [0058] Finally, the processor produces 103 the saliency map based on the image data and the distorted image data. That is, the processor calculates the saliency map by performing a mathematical operation using the image data and the distorted image data as inputs and the saliency map as an output. In that sense, the calculation may be pixel- wise, so that one output pixel of the saliency map is a result of a calculation using one pixel of the image data and one pixel of the distorted image data, where all pixels are taken from the same image location. In example, the saliency map may be a comparison between the image data and the distorted image data. In this sense, the pixel value at each respective location of the image data and the distorted image data may be compared. The comparison may be the absolute value of the difference in pixel intensity values at each respective location. In another example, the saliency map may be a subtraction of the distorted image data from the image data. In this sense, the pixel value of each respective location in the image data may be subtracted from the pixel value of the distorted image. Therefore, in this example, the subtraction between the pixel value at each respective location would result in a value that possesses both an absolute value of the subtraction, but also a sign of the subtraction. [0059] The difference between the image data and the distorted image data contrasts the difference between the first and second classification, which is seen visually in the saliency map. Further, the sign of the subtraction between the image data and the distorted image data, as well as, the absolute value of the subtraction may be indicative of physical aspects captured in the image data. For example, if the image data represents a MR image of a patient’s brain and the image data is subtracted from the image data, areas where subtraction is negative may be indicative of tissue atrophy (decrease in size of a body part, cell, organ, or other tissue). Further, areas where the subtraction is positive may be indicative of an increase in the tissue, which may represent tissue swelling. Even further, the absolute value of the subtraction may be indicative of the amount of tissue atrophy or increase. The sign of the subtraction, as well as the absolute value, may be colour-coded to make visually clearer to the user and to make interpretation of the subtraction easier. [0060] The saliency map indicates salient features of the image data in relation to the classifier. Salient features may also be referred to as important features, which the classifier deems important or significant when determining the classification of image data. For example, if the classifier is used in the medical profession for the diagnosis of a disease, the salient features may indicate symptoms associated with the diagnosis of the disease. These symptoms may include swelling of the brain, for example. As such, the disclosed method may be used to better understand or determine symptoms with a diseases, if symptoms of a disease are understood. The saliency map may also be superimposed on the image data to highlight the salient features within the original image data. [0061] Fig.2 illustrate an example of a spatial-sequence 200 of two-dimensional cross-sectional images. Both the image data and the distorted image data may comprise such a spatial sequence. It is this spatial sequence of two-dimensional cross sectional images that forms the input into the classifier. The classifier then determines the first or second classification based on the spatial sequence of two-dimensional cross sectional images used as input. The spatial sequence 200 comprises a spatial dimension 201 and a pixel value 202 of one of the two-dimensional cross-sectional images. In Fig.2, the spatial dimension 201 is indicative of the overall spatial sequence 200 of two- dimensional cross-sectional images. Here, the spatial dimension 201 is represented by the spatial variable z. The remaining two-dimensions are represented by the spatial variables x and y. The pixel value 202 is indicative of a Hounsfield value in a slice of a MRI or CT scan. [0062] The spatial sequence 200 of the two-dimensional cross-sectional images represent a scan volume comprising a plurality of voxels, each of the plurality of voxels has a voxel value. The voxel values are based on the pixel values from the two- dimensional cross-sectional images in the spatial sequence. In this sense, a voxel can be considered as a 3D pixel and the voxel value is based on, and may be identical to, a pixel value. In one example, a voxel is defined by the overlapping pixel between two adjacent cross-sectional images and the corresponding voxel value is based on one of the values of the overlapping pixels. In another example, the voxel value is based on both pixel values of the overlapping pixels, such as an average of the two pixel values. In a further example, a voxel is defined by the overlapping pixel between three or more adjust cross-sectional images and the corresponding voxel value is based on any one of one the values of the overlapping pixels, such as an average of the pixel values. In this way, each voxel value of the plurality of voxels corresponds to a Hounsfield value. [0063] If the image data is indicative MRI scan, then the image data may be created using a MRI scanning machine. MRI scanning machines consists of a table that slides into a cylinder. Inside the cylinder is a magnet that creates a strong magnetic field within the cylinder. The strong magnetic field is used to align the protons in the water molecules in a patient’s body, as the protons possess a magnetic moment due to their intrinsic spin. Once the protons align with the strong magnetic field, a radio frequency signal is sent from the MRI scanning machine, which excites the protons and causes them to align against the strong magnetic field. Once the radio frequency signal is turned off, the protons relax and re-align with the strong magnetic field, thereby emitting electromagnetic radiation, which is detected by the MRI scanning machine. A patient lies on a motorised table that moves horizontally, in increments, through the cylinder and captures a single MR image at that increment, which corresponds to a two- dimensional cross-sectional image (virtual “slice”) of the patient’s body. Each two- dimensional cross-sectional image corresponds to a different region of the patient’s body as they move horizontally, in increments, through the cylinder. The collection of two-dimensional cross-sectional images also has an associated spatial sequence, in which the cross-sectional images are ordered along. The horizontal direction that the patient moves within the cylinder is indicative of a dimension in space along the spatial sequence of the two-dimensional cross-sectional images. The collection of two- dimensional cross-sectional images provides a three-dimensional representation of a patient’s internals, which aids physicians in diagnosing internal diseases and conditions. [0064] The two-dimensional cross-sectional image consists of black and white regions, as well as shades of grey. These regions correspond to the rate at which excited protons return to their equilibrium state after being exposed to radio frequency (equilibrium state corresponds to the spin state which aligns with the strong external magnetic field), as well as the amount of energy released, changes depending on the environment and the chemical nature of the molecules. These factors can be differently weighted in the MR image, to give more contrast to certain tissues in the image. For example, in a T1-weighted (T1w) image, regions of high signal representing white regions may correspond to protein-rich whereas in T2-weighted (T2w) images, these regions may correspond to water rich fluids in the patient’s body. [0065] If the image data is indicative of a CT scan, then the image data may be created using a CT scanning machine. CT scanning machines use a rotating X-ray tube and a row of detectors placed in a gantry of the machine to measure X-ray attenuations by different tissues inside the body. X-ray attenuations refers to when X-ray photons, with wavelengths ranging from 10 picometers to 10 nanometers, are absorbed by the tissue when the X-ray photons pass through the body of a patient. A patient lies on a motorised table that moves horizontally, in increments, through the rotating X-ray tube and captures multiple X-ray measurements from different angles at each increment. The multiple X-ray measurements are then processed on a computer using reconstruction algorithms to produce tomographic two-dimensional cross-sectional images (virtual “slices”) of a body. Each two-dimensional cross-sectional image corresponds to a different region of the patient’s body as they move horizontally, in increments, through the rotating X-ray tube. Similarly to MR images, the collection of two-dimensional cross-sectional images also has an associated spatial sequence, in which the cross-sectional images are ordered along. The horizontal direction that the patient moves within the cylinder is indicative of a dimension in space along the spatial sequence of the two-dimensional cross-sectional images. Similarly to MR images, the collection of two-dimensional cross-sectional images provides a three-dimensional representation of a patient’s internals, which aids physicians in diagnosing internal diseases and conditions. [0066] The two-dimensional cross-sectional (X-ray) images consists of black and white regions, as well as shades of grey. These regions correspond to the areas of the body that have attenuated the X-rays. For example, structures such as bones readily absorb X-rays and, thus, produce high contrast on the X-ray detector. As a result, bony structures appear whiter than other tissues in the X-ray image. Conversely, X-rays travel more easily through less radiologically dense tissues such as fat and muscle, as well as through air-filled cavities such as the lungs. These structures are displayed in shades of grey in the X-ray image. where little to no attenuation has occurred will show up as black regions in the X-ray image. [0067] Fig.3 illustrates a system 300 for creating a saliency map. The system 300 comprises a computer 301 which itself comprises a processor 302 connected to a program memory 303 and a data memory 304. The program memory 303 is a non- transitory computer readable medium, such as a hard drive, a solid state disk or CD- ROM. Software, that is, an executable program stored on program memory 303 causes the processor 302 to perform the method in Fig.1, that is, processor 302 applies the classifier to the image data to determine the first classification, distorts the image data to produce distorted image data and produces the saliency map. The data memory 304 may store the classifier. For example, if the classifier is a trained neural network, the data memory 304 may store the classifier in the form of network parameters that have been optimised for determining a classification of the image data. The data memory 304 may store image data indicative of an MRI scan or CT scan as a single-file DICOM (Digital Imaging and Communications in Medicine), multi-file DICOM format or a similar/equivalent image format, for example. The image data may also be indicative of a two-dimensional image, which may be stored on data memory 304 as Joint Photographic Experts Group (JPEG) format, RAW image format or a similar/equivalent image format, for example. [0068] The processor 302 may receive data through all these interfaces, which includes memory access of volatile memory, such as cache or RAM, or non-volatile memory, such as an optical disk drive, hard disk drive, storage server or cloud storage. The computer system 300 may further be implemented within a cloud computing environment, such as a managed group of interconnected servers hosting a dynamic number of virtual machines. [0069] It is to be understood that any receiving step may be preceded by the processor 302 determining or computing the data that is later received. For example, processor 302 may then store the saliency map on data memory 304, such as on RAM or a processor register. Processor 302 then requests the data from the data memory 304, such as by providing a read signal with a memory address. The data memory 304 provides the data as a voltage signal on a physical bit line and the processor 302 receives the image data as the input image via a memory interface. [0070] Fig.1 is to be understood as a blueprint for the software program and may be implemented step-by-step, such that each step in Fig.1 is represented by a function in a programming language, such as C++ or Java. The resulting source code is then compiled and stored as computer executable instructions on program memory 303. [0071] It is noted that for most humans performing the method 100 manually, that is, without the help of a computer, would be practically impossible. Therefore, the use of a computer is part of the substance of the invention and allows performing the necessary calculations that would otherwise not be possible due to the large amount of data and the large number of calculations that are involved. [0072] System 300 may also comprise an imaging machine 305, an internet server 306 and a monitor 307, which is in communication with the computer 301 via an input/output port 308. In an example, imaging machine 305 may be a computed tomography (CT) scanning machine or an MR imaging (MRI) machine, for which the image data would be indicative of a CT scan or a MRI scan. In the example where imaging machine 305 is a CT scanning machine, the imaging machine 305 scans a patient by incrementally capturing two-dimensional cross-sectional images of the patient. Each of the two-dimensional cross-sectional images are indicative of an X-ray image, which correspond to a different region of the patient’s body as they move horizontally, in increments, through the CT scanning machine. The collection of cross- sectional images represent image data that is indicative of a CT scan and represents a 3D scan of the patient. Once CT scanning machine creates the image data of the patient, processor 302 receives the image data via the input/output port 308. A similar process occurs in the example, if the imagining machine 305 were a MRI machine. Processor 302 then performs method 100 on the image data received from the imaging machine 305 and creates a saliency map from the image data. Additionally, the image received from the imaging machine 305 be stored on the data memory 304 before and after method 100 is performed by processor 302. [0073] In another example, the internet server 306 may store image data and/or the distorted image data. Additionally, image data and/or the distorted image data may be communicated to the internet server 306 and the internet server may then communicate the image data to computer 301. The image data may be received from a source external to system 300, such as another system located remotely to system 300. The remote system may be located within the same facility as system 300 or may be completely remote from system 300. The image data may be stored in the internet server 306 as a single-file DICOM or multi-file DICOM format. The internet server 306 may be in communication with the imaging machine 305, or may be in communication, over a communications network. Similarly, processor 302 receives the image data via the input/output port 308 and performs method 100 on the image data. After processor 302 performs method 100 on the image data, the image data and/or the saliency map may be communicated to the internet server 306. The internet server 306 may then communicate the image data or the saliency map to an external system or store them. [0074] System 300 may also comprise a monitor 307, which may be configured to display the image data from the imaging machine 305 and from the internet server 306. In this way, processor 302 communicates the image data from the imaging machine 305 and from the internet server 306 to the monitor via the input/output port 308. Monitor 307 may be further configured to display the saliency map and/or the distorted image data after processor 302 performs method 100 of the image data. In this way, processor 302 communicates the saliency map and/or the distorted image data to the monitor via the input/output port 308. [0075] Software may provide a user interface that can be presented to the user on a monitor 307. The user interface is configured to accept input from the user, via a touch screen or a device attached to monitor 307 such as a keyboard or computer mouse. These device may also include a touchpad, an eternally connected touchscreen, a joystick, a button, and a dial. In an monitor 307 may display multiple instances of image data, and the user may choose one of the multiple instances of image data for method 100 to be performed on, by processor 302. The user may choose one of the multiple instances of image data by interacting with the touch screen or inputting the selection with a keyboard or computer mouse. The user input may also be provided to the input/output port 308 by the monitor 307. Processor 302 may then perform method 100 on the selected image data and communicate the image data and/or the saliency map to monitor 307. [0076] The processor 302 may receive data, such as image data and/or the distorted image data, from data memory 304 as well as from the communications port 308. In one example, the processor 302 receives image data from the imaging machine 305 and internet server 306 via communications port 308, such as by using a Wi-Fi network according to IEEE 802.11. The Wi-Fi network may be a decentralised ad-hoc network, such that no dedicated management infrastructure, such as a router, is required or a centralised network with a router or access point managing the network. [0077] Although communications port 308 is shown as single entity, it is to be understood that any kind of data port may be used to receive data, such as a network connection, a memory interface, a pin of the chip package of processor 302, or logical ports, such as IP sockets or parameters of functions stored on program memory 303 and executed by processor 302. Computer 301 may be configured with multiple input/output for each of the imaging machine 305, internet server 306 and monitor 307. The parameters of functions may be stored on data memory 304 and may be handled by-value or by-reference, that is, as a pointer, in the source code. [0078] Fig.3 is one example of a configuration of system 300. However, system 300 is not strictly limited to this configuration and this may be one possible embodiment of system 300. In another embodiment, computer 301 may not be in direct communication with imaging machine 305 or monitor 307 and may only in communication with the internet server 306. Method [0079] The terms distorted, warped and counterfactual are also used interchangeability to describe the image data that is distorted and gives the opposite classification to the original image data, thereby corresponding to the second classification. [0080] The disclosed method is described with reference and applications to medical imaging. In this sense, image data may comprise images of an anatomy of a patient. An example implementation of the disclosed method is presented here for creating a saliency map for an MR image and the classifier determines the presence of Alzheimer’s disease. However, this disclosed method is not limited to MR images, nor it is limited to classifiers of Alzheimer’s disease or neurological diseases. Further, the disclosed method is not limited to creating saliency maps for the purposes of medical imaging, but may be applicable for any classifier that classifies image data. This is only one example to demonstrate the disclosed method for creating a saliency map for the classification of image data. [0081] The disclosed method may also be performed on two-dimensional (2D) image data or three-dimensional (3D) image data. In this case where the image data is 3D, the image data may comprise images that form a spatial sequence of two-dimensional cross sectional images. In this example, the image data may be indicative of a CT scan or an MRI scan. In other words, the image data comprises magnetic resonance image data or computed tomography image data. In this case, the spatial sequence of two-dimensional cross sectional images would form an input to the classifier to determine the first or second classification. [0082] The classifier may be a machine learning model, such as any variation of neural network including convolutional or deep, decision tree or support vector machine (SVM). In an embodiment, the classifier is a trained machine learning model. These are only a few examples of possible machine learning models that the classifier may be, however the classifier is not limited to the models mentioned here. The classifier may also comprise a of machine learning models, in which the output of one machine learning model may be input for a subsequent machine learning model. In a preferred embodiment, the classifier has exactly two possible output classifications. In this sense, the classifier may be a binary classifier. [0083] In an example, the classifier is a convolutional neural network (CNN), which is ideal for applications involving images and image data. Further, in this example, the classifier is a three-dimensional (3D) CNN, in which the convolutional filter moves in three dimensions, as opposed to two dimensions in a standard CNN. A 3D CNN is ideal for applications involving image data that comprise multiple images, such as where the image data is indicative of MRI or CT scan. In this example, the convolutional filter moves along the horizontal and vertical directions of each individual images, as well as along the direction that is indicative of the spatial sequence, which may also be referred to as the dimension defined in space along the spatial sequence. [0084] Mathematically, a convolution is an integration function that expresses the amount of overlap of one function g as it is shifted over another function f. Intuitively, a convolution acts as a blender that mixes one function with another to give reduced data space while preserving the information. In terms of Neural Networks and Deep Learning, convolutions are filter (matrix / vectors) with learnable parameters that are used to extract low-dimensional features from an input data. They have the property to preserve the spatial or positional relationships between input data points. CNNs exploit the spatially-local correlation by enforcing a local connectivity pattern between neurons of adjacent layers. In this example, the CNN architecture may comprise an optimal number of convolutional layers, filter size and stride length. These factors may be optimal for determining the first and second classification. [0085] Intuitively, convolution is the step of applying the concept of sliding window (a filter with learnable weights) over the input and producing a weighted sum (of weights and input) as the output. The weighted sum is the feature space which is used as the input for the next layers. [0086] In this example, the CNN may also comprise a batch normalization layer and ReLU activation function. Further, the output of the last convolution layer may be subjected to Global Average Pooling (GAP). Even further, a fully connected layer with Sigmoid activation may also be used on the output layer. However, other activation functions may be used. These activation functions include, but are not limited to, a binary step function, a tanh function, a ReLU function or a softmax function. [0087] Fig.4 illustrates a 3D convolution applied to the data structure from Fig.2.3D convolutions apply a 3D convolutional filter to the dataset and the filter moves in 3- directions (x, y, z) to calculate the low level feature representations. Their output shape is a 3 dimensional volume space such as cube or cuboid. [0088] In Fig.4, a filter 401, which is a 3 x 3 x 3 filter in this example, calculates a weighted sum of voxel values. There are 9 weights, which are the coefficients of the filter. So the filter starts at an initial position 402 in the data structure 200, multiplies each voxel value with the respective filter coefficient and adds the results. Finally, the filter 401 stores the resulting number in an output pixel 403. In this sense, the output voxel values of each of the one or more convolutional layers comprises a weighted sum of input voxel values. The weights in the weighted sum correspond to the coefficients of the filter. Then, the filter moves by one voxel along one direction in the data structure 200 and repeats the calculation for the next voxel of the output image. For example, if the stride is 2, the filter will move by two voxels along one direction. That direction may be an x-dimension 404 or a y-dimension 405. Importantly, the filter 401 may also move in space along the spatial sequence 406. That is, the cubic 3 x 3 x 3 filter 401 adds voxel values along the spectral dimension 406 and also moves along the spectral dimension 406. Therefore, this operation is referred to as a convolution along a dimension defined in space along the spatial sequence. It is noted that there are typically more than one filters 401, such as two or more, which create respective 3D output data structures (only one is shown in Fig.4). [0089] A single convolutional layer multiple filters, where each filter corresponds to a different feature. The corresponding feature of each filter is determined during training of the CNN. More specifically, each filter may contain weights, which are adjusted during a training process through a backpropagation process, for example. The filters are used to quantitatively determine the contribution of a particular feature on the output of the CNN. For a convolutional layer that occurs at the start of a CNN architecture, the filters of this convolutional filter may correspond to simple features. However, convolutional layers that occur later in the CNN architecture may exhibit more complex or abstract features. As an example, the complex or abstract features may be a combination of the simple features from the previous convolutional layer. [0090] Fig.5 illustrates an output layer 500 of the classifier, in the example where the classifier is a CNN, which is fully connected. The output layer has inputs 501 and one output 502. The inputs 501 are the results from the previous convolutional layers. That is, the convolutional layer, such as in Fig.4 is connected to the output layer 500 to provide the results to the output layer 500. However, the convolutional layer need not be connected directly to the output layer 500 in the sense that there may be further non- convolutional layer between the convolutional layer and output layer 500. For example, there may be a max-pool layer or a GAP layer between the convolutional layer and the output layer 500. [0091] The output layer has and output 502 and provides an output value that is indicative of the first or second classification. It is noted that Fig.5 shows a simplified example of only four inputs 502, noting that in most practical applications, the number of inputs is significantly higher. [0092] The processor 302 may also train the classifier by calculating the error between the prediction of given image data and the actual classification of the image data, then minimising this error. If the classifier is a CNN, for example, the error between the prediction (or the predicted classification) and the actual classification can be minimised by adjusting the weights associated with the weight sum of the convolutional filter. The processor 302 these weights through the process of back propagation and gradient descent to gradually improve the weights of the CNN and reduce the error. [0093] Training data is used to train the classifier, which comprises label data, such as image data that has been labelled with a classification, where the label for each item of training data may be manually provided by a user. For example, if the classifier classified image data as having Alzheimer’s disease or not having Alzheimer’s disease, then the training data may be image data, such as MRI scans or CT scans, that have been diagnosed by a medical or health professional, as having or not having Alzheimer’s disease. Training data may also or otherwise be obtained from a database containing already labelled data, for which manually label the data is not needed. Training the classifier in this way is known as “supervised learning”, in the sense that, the prediction of the classifier is validated and the weights of the classifier are adjusted accordingly. It should be noted that “supervised learning” is not the only method that can be used to train the classifier. For example, “semi-supervised learning” can be used in situations where only a relatively small amount of training data can be obtained. [0094] Fig.6a illustrates an overview of the disclosed method. Given a classifier, which may be a trained machine learning model such as a Deep Neural Network (DNN), an input image x and its prediction ^p , the disclosed method seeks to optimize a set of FFD grid-points ^ such that when the distorted image x ^ is fed to the classifier, the prediction is swapped, i.e., applying the DNN to the distorted image x ^{^} gives the opposite prediction to p (which can be represented as p ). In this disclosure, the classifier is a DNN, for explanatory purposes, so the classifier will often be referred to as a DNN in this disclosure. Here, the prediction p may be referred to as the first classification and the opposite prediction p may be referred to as the second classification. Given a DNN _{f ^} , processor 302 smoothly deforms an input image x to produce a counterfactual example x ^ . The distorted image is counterfactual, in the sense that, the prediction may by the classifier is opposite to the true classification obtained by applying the classifier to the original image data. [0095] Fig.6b illustrates the original x , distorted image x ^ and the corresponding saliency map x ^{^} ^ x . Producing the saliency map is based on the image data x and the distorted image data x ^ and comprises calculating a difference between the image data and the distorted image data. The saliency map can then be superimposed onto the original image x to highlight the salient features of the image that are important to the DNN when making the prediction. [0096] More specifically, consider a trained DNN f ^{H ^ W ^ D} ^ : ^ ^ ^ that uses a 3D image as input and outputs a real number, where the corresponds to the classification of the DNN. Now consider the DNN f _^ as well as an input imagex ^ ^ ^{H ^ W ^ D} for which one wants to provide a saliency map for the decision p . Given ^ , x ) , processor 302 optimizes a smooth deformation (distortion) ^ _^ ( ^ ) to _{produce a counterfactual prediction f ^ (x} ^{^} _{) ^ p from the deformed image} ^{( x )} _{. In} other words, creating the saliency map method comprises determining of the image data. The distortion ^ _^ ( ^ ) can be considered as an operator that acts on image data x . The smooth ^ _^ ( ^ ) can be adjusted by changing the set of grid points ^ . The set of grid points provides a counterfactual or opposite predict may be denoted as ^ ^{^} . In other words, determining a distortion of the image data comprises adjusting the distortion to optimise an output of the classifier applied to the distorted image data. The output of the classifier can be considered “optimised” when the second classification is opposite to the first classification. In other examples, “optimised” refers to the smallest distortion that is necessary to achieve an opposite classification. [0097] The distortion of the image data is smooth and invertible, while the distorted image produced as a result of the distortion, is anatomically plausible, when the disclosed method is applied in the field of medical imaging. Anatomical plausibility of the distorted image may be viewed as the inability to anatomically distinguish the distorted image from a captured image of anatomy. In the disclosed method, processor 302 creates a distorted image that is similar to the original image. The images are similar, in the sense that they would be from visual inspection. As the original image would be anatomically plausible, as it would be indicative of the anatomy of a patient, the disclosed method produces an anatomical plausible distorted image. Parameterisation of the deformation. [0098] Here, the equations and following method apply to 3D image data. However, the disclosed method is not limited to 3D image data. The disclosed method may be similarly applied to 2D image data through simple modification of the equations and following method, for example, reducing the dimensionality of the below equations. To produce smooth deformations that are invertible and anatomically plausible, the FFD is parameterised by a grid of points ^ ^ ^ ^{N ^ N ^ N ^ 3} . In some examples, N is a factor of the dimension of the image data. For example, the image data may have a dimension of H x W x D = 256 x 256 x 256 = 256 ³ , and N is a factor of 256, such as 32. ^{[0099] It is noted the space complexity of finding an adversarial deformation} _{^ ^ ( ^ )} may depend on the size of N . However, N should set to a large enough value to avoid deformations that are too coarse. The final x 3 in the dimensionality of ^ stems from the spatial nature of the deformation. [0100] More generally, ^ _^ ( ^ ) is the re-sampled digital image determined by the deformation vector ³ ^ ^ ³ which is defined as, v _^ ^{( x} ₁ ^{, x} ₂ ^{, x} ₃ ^{) =Id} _{^ 3} ^{( x} ₁ ^{, x} ₂ ^{, x} ₃ ^{) ^ u} _^ ^{( x} ₁ ^{, x} ₂ ^{, x} ₃ ^), ₍₁₎ [0101] with Id the identity map, and with _{u ^} defined as, 4 ^{^} n ^ _^ ^{^ ^ x} ^ i ^{^ ^ x} i ^{^ ^ ^} ^ ^{^ (2)} [0102] where ^ ^ ^ denotes the floor ^ ^{^} _i = ^ ^xi ^ 1 ^ _{, and N} ^{^} _^ ^th i _{denotes the ^i} cubic B-spline basis functions, ^ ^{t = 1} 6 ^{1 ^ t 3 ^} ₂ ^{t = 1} 6 ^{3 3 2} 1 ^{^ ^ ^ ^ ^ ^} _^ ^{t ^ 6 t ^ 4} _^ t ^. [0103] It is noted a ^ ³ , while ^ represents a global deformation. In essence, both ^ and ^ represent the same concept, however they are two different objects. [0104] In essence, the distortion is parameterised by a set of grid points that span the image data. If the image data were 3D, then the image data would comprise a plurality of voxels, each with a voxel value, and distorting the image data with the free-form deformation comprises applying a deformation vector field to voxels of the image data. Similarly, if the image data were 2D, then the image data would comprise a plurality of pixel, each with a pixel value, and distorting the image data with the free-form deformation comprises applying a deformation vector field to pixels of the image data. The pixel or voxels value may be indicative of a Hounsfield value used in the medical imaging field. Optimisation algorithm. [0105] Given the aforementioned FFDs, the optimisation problem is described as: 3 ^ ^{^} = argmin _{^ ^} x _^ ^ _^ x ^{^} ^ ^ ^{^} ), ^ _^ ⁽ x ⁾ _^ ^ ^{( f} ^ ⁾ [0106] where ^ (f _^ ) is the decision for the DNN f _^ that one wants to ^{reach. Without loss of generality, suppose that} _{f : ^} ^{H ^ W ^ D} ^ _{^ [ ^ 1,1]} ^{is a binary} classifier, the resulting decision boundary reads as ^(f _^ ) = ^x : f _^ ( x ) = 0 ^ . The sign of the first classification is used in the above so that the optimisation may work regardless of whether is -1 or 1. In essence, this optimisation seeks to find the set of grid points ^ which causes the product sign(f _^ (x )) f _^ ( ^ _^ ( x )) to change sign. In essence, optimising an output of the classifier applied to the distorted image data further comprises repeating the optimisation until the output of the classifier indicates that the second classification is opposite to the first classification. [0107] ^ _q in the optimisation problem represent weighting factors, which weigh the contribution of the penalty terms ^ _q . In other words, ^ _q represent parameters (or hyperparameters), which can be tuned to quantify of each penalty term ^ _q . Processor 302 may determine or tune these parameters by conducting a grid search or a parameter sweep to find an optimal set of ^ _q that enable the optimal overall penalty to the above optimisation problem. For example, each ^ _q may be a value between 0 and 1, and processor 302 may calculate the above optimisation problem for each value of ^ _q between 0 and 1 at a predetermined step size. Processor 302 may then parameter set, which results in the optimal overall penalty to the above optimisation problem and return that parameter set as the solution to the grid search. [0108] In some examples, ^ ₁ and ^ ₃ are parameters that allow for an anatomically plausible image, in the sense that these parameters tune the contribution of penalty terms that enable the deformation to be anatomically plausible. As such, ^ ₁ may prevent the warp from creating folding that results in artefactual intensity change in the image (non-clinically plausible). Further, ^ ₃ may enable the user to prescribe the deformation to a specific anatomic area. ^ ₂ may be a term, as several deformations fulfilling the requirements may be found. Hence, ^ ₂ may be tuned to find the deformation with minimal energy. In other words, ^ ₂ may be used to find the deformation with the minimum amount of effort. [0109] In light of this, in some examples, ^{^} ₁ and ^{^} ₃ may be the most important parameters, as ^ ₁ enables the deformation to be invertible and ^ ₃ restraints the transformation to a specified region. As such, these parameters allow the counterfactual image to be anatomically plausible, while focusing only the anatomy area of interest. Thus, in this example, ^ ₁ and ^ ₃ would be greater than ^ ₂ , which results in the second penalty term contributing less to the overall penalty in the optimisation problem. As such, in the example where the classifier classifies image data as having (or not having) Alzheimer's disease, these parameters may be tuned to enable the deformation to change only cortical structures and not affect the patient's skull. [0110] For explanatory purposes of the optimisation problem, consider the case where the first classification is given by f _^ (x ) ^ ^ 1. Before applying the optimisation algorithm, the distortion ^ _^ would be equal to the identity transformation i.e., ^ _^ (x ) ^ x and therefore, the second classification would be equal to the first classification initially i.e., f _^ ( ^ _^ (x )) ^ ^ 1. Thus, the product sign(f _^ (x )) f _^ ( ^ _^ ( x )) would be equal to 1 initially. However, the optimisation algorithm would stop when f _^ ( ^ _^ ( x )) reaches the opposite classification. Therefore, the product sign(f _^ (x )) f _^ ( ^ _^ ( x )) would be equal to -1 when the distorted image reaches the opposite classification to the original image and hence, the optimisation algorithm stops when the product sign(f _^ (x )) f _^ ( ^ _^ ( x )) changes sign. Similar is true if the first classification is 1, rather than -1 in the above example. [0111] Processor 302 determines a solution to the optimisation problem using a gradient descent approach, as described by Algorithm 1. This can be seen as optimising a cost function and calculating a update the grid points in the direction of the gradient, wherein the cost function is a function of the distortion. The cost function is indicative of the argument of the argmin function. In Algorithm 1, processor 302 repeats the optimisation until the output of the classifier indicates that the second classification is opposite to the first classification. Algorithm 1 Input Imagex , model f _^ Output Image x ^ Initalizex ₀ ^x, i ^ 0, ^ ₀ ^ 0 _{^ N ^ N ^ N ^ 3} while f _^ ^ ^{^ sign( f (x )) f r} i _{^ ^} ^{( x} _i ^{) ^} ^ _q =1 ^{^} _q ^{^} _q ^{( ^} _i ⁾ step x _{i ^ 1} ^ ^ _{^i ^ 1} ( x _i ) Image update i ^ i ^ 1 end while return x _i [0112] As can be seen from Algorithm 1, determining the distortion ^ _^ ( ^ ) comprises iteratively updating the set of grid points until the output indicates that the second classification is opposite to the first classification. Algorithm 1 proceeds until the ^{distorted image} _x ^{^} _{^ f ^ ( x )} ^{crosses the decision boundary of the classifier. In this sense,} the distorted image data is indicative of image data on a decision boundary of the classifier. Penalties [0113] It is assumed that disease can be modelled by smooth and invertible mappings and the distortion from the healthy control to the disease group can be generated using a sufficiently refined FFD. To ensure these assumption are upheld, constraints are enforced on the distortion ^ _^ , which is implemented by adding penalty terms to the cost function. In this the distortion _{^ ^} comprises at least one penalty term to enforce a constraint on the distortion. [0114] The first restriction considered is the local invertibility of the distortion provided by the FFD. The quadratic penalisation is given by: 1 3 N ^ ¹ ₍ ^{^} ₎ ^{^} ^ ^ p ₍ ^{^} i _^ ^1, j ^, k ^, d ^{^ ^} i j k d _; ^ d , x _, ^ d , x ₎ N ₃ ^{, , , 1 2} ^ 1 ^ 1 3 ), [0115] with p a quadratic penalty function defined as, ^ ^{1 ^} ₍ t ^{^} ^ 1 ₎ 2 _if t _< ^ 1 [0116] Local invertibility allows the user to go backward, i.e. given the counterfactual image, to find a deformation of the voxels that result in recreating the original image. This is useful when creating counterfactual images as it prevents the creation of foldings, which prevents the counterfactual image from being anatomically plausible. A folding, in this context, can be considered to be when the image overlaps itself (i.e., pixel or voxel values add onto one another) during a deformation of the image, rather than areas of the image moving during the deformation. [0117] The second constraint enforces the sparsity of the distortion, [0118] The sparsity of of ^ ₂ among all of the anatomically plausible counterfactual (adversarial) deformations. As previously discussed, ^ ₂ may be used to find the deformation with the minimal amount of energy. [0119] Finally, the support of ^u is restrained to ^ ^ ^ ³ , ^ ³ ₍ ^{^} _{) =} ^{1 ^} ₍₁ ^{^ ^} ₎ ^{^ ^ ^} _{, (4)} N ₃ ^{^ 1} [0120] where ^ _^ a discrete binary mask down-sampled to the grid points’ resolution. In this sense, this penalty term is indicative of a subset of the image data, which restricts the distortion to only target the anatomy area of particular interest. This can aid the classifier in making its classification, as the classifier may only be trained with image data associated with a particular anatomic area. [0121] For illustration and explanatory purposes, ^ is restricted to the brain in this disclosure, as the disclosed method is presented using a neural network trained to classify neurodegenerative diseases from brain MRI. Implementation details. [0122] As previously indicated, the disclosed method is model agnostic and is not only restricted to convolutional architectures. In an example, the disclosed method efficiently implements the FFD by leveraging strided transposed convolutions resulting in a small memory footprint of 4.92 GiB and a fast runtime of 0.256 ^ 0.008s per gradient step. Algorithm 1 converges in 32.9 ^ 13.4 iterations ( ^ 8.45 s) for a 256 ^ 256 ^ 256 image. Experiments [0123] To demonstrate the effectiveness of the disclosed method, the disclosed method is evaluated using a DNN, which, given a three-dimensional T1w structural MR image, predicts if the patient is healthy or has Alzheimer’s disease. The architecture is inspired by a lightweight convolutional neural network (CNN). This model is trained using early stopping on an augmented version of the Alzheimer’s disease Neuroimaging Initiative (ADNI) dataset with a train/validation/test split of 776/342/520 (the ADNI dataset is available at: https://adni.loni.usc.edu/). It is noted that different splits, corresponding to different IDs and time points from the ADNI dataset, are used to evaluate the model. [0124] It should be noted that training the classifier is not a necessary step in the disclosed method. The disclosed method is applicable to any machine learning model that has been trained for purposes unrelated to the disclosed method. However, for the purposes for evaluating the disclosed method in comparison to other saliency map methods, it was preferred to train the classifier. Details of the DNN [0125] Fig.7 illustrates the architecture of the DNN used to evaluate the disclosed method. Training the DNN comprised perform a grid search to find the best set of hyperparameters using an Adam and/or SGD optimiser with the following learning rates: {5e−1, 1e−1, 5e−2, 1e−2, 5e−3, 1e−3, 5e−4, 1e−4, 5e−5, 1e−5, 5e−6, 1e − 6}. [0126] The mini-batch size was set to 8. The best model was selected using early stopping; the resulting hyperparameters were the Adam optimizer and a learning rate of 5e − 5. During training, the dataset is augmented using TorchIO (https://torchio.readthedocs.io/) which is an open-source Python library for efficient loading, preprocessing, augmentation and patch-based sampling of 3D medical images in deep learning, following the design of PyTorch. Here, the following augmentations were made use of: ^ RandomAffine; ^ RandomBiasField; and ^ RandomNoise which are defined from line 64 to 66 of Code/code dnn/data.py. On the test set, the selected model has an accuracy of 83.53% a sensitivity of 83.43% and a specificity of 80.70%. Healthy Control Alzheimer’s Disease S N Mean A e (± % S N Mean A e % 7 7 2 5 of subjects and data po nts respect ve y. Resources [0127] The neural networks were trained using a cluster composed of 112 nodes (each node is equipped with 4 P100 GPUs (16 GB) 2×14 CPUs and 256GB of RAM). Each model was trained using 4 CPUs and 16GB of RAM and 1 GPU. Robustness [0128] Fig.8 illustrates robustness assessment of the disclosed method for Gaussian noise with a standard deviation of 0.1 injected on the input image. The robustness of disclosed method is evaluated towards random perturbation of the input image (Gaussian noise). For a standard deviation of σ = .1 (10%, images are normalized within [0, 1]), the resulting dense deformation field norm is displayed for six different perturbed images. A new visualisation technique provided disclosed method [0129] The counterfactual images x ^{^} provided by Algorithm 1 are very close to the initial images x . Indeed, the global deformation of the FFD defined in Equation (1) is set to be the identity. Moreover, the regularisation term in Equation (3) enforces the deformation to be sparse. As a result, the difference of x ^ x ^{^} for a given patient can be used to display voxel intensity variations between the original and the delusive image accurately. The disclosed method appears to modify the size of the ventricles and the thickness of the temporal cortices. This pattern tends to be reproducible across the studied cohort; a global visualization is proposed based on registration to access the comprehensiveness of the observed behaviour. _{[0130] Fig. 9 illustrates the saliency map x} ^{^} _{^ x overlayed on the on the initial} patient’s MR Image. This image depicts the slight change of the initial image that results in the model changing its prediction. It is noted that some areas have a positive value (increase in brain mass) while other areas have a negative value (reduction/atrophy) in brain mass. In this example, the image may be classified as having Alzheimer and it should be noted that the areas with both the positive and the negative values have contributed positively to (i.e. support) the Alzheimer classification. That is, it is not the case that the negative values indicate no Alzheimer and the positive values indicate Alzheimer. Instead, both areas of negative and positive values are salient in the sense that they lead towards an Alzheimer classification. Quantitative analysis of the disclosed method’s saliency map and comparison to pre- existing methods [0131] The saliency maps derived from the disclosed method are now compared against three techniques for post-hoc interpretability (Implementation recovered from TorchRay https://facebookresearch.github.io/TorchRay/). [0132] Fig.10 illustrates saliency maps for each method of an axial and a sagittal view given a single prediction. Percentage of total saliency map averaged for all of the scans present in the experimental split correctly classified as Normal Control (NC). In Fig.10, saliency maps are depicted for each of the methods when considering images that are correctly classified as a NC. Saliency maps for DeconvNet and Gradient are very scattered and diffused, thus challenging to interpret. On the contrary, maps issued from GradCAM are very evasive, highlighting a large zone that leaves room for interpretation bias. Finally, these methods are compared by computing the averaged percentage of energy present in the patient’s cortical sub-structures using FreeSurfer’s atlas. The saliency map derived from the disclosed method indicates that most of the deformations are happening in the temporal lobes, hippocampal area and ventricles, which are discriminating regions for diagnosing of Alzheimer’s disease (AD). Robustness to localised deformations [0133] The robustness of a prediction to extremal deformations carried out on individual anatomical substructures is now tested. Using the segmentation masks of FreeSurfer, Algorithm 1 is modified to restrict the deformation sequentially to every parcellation from this atlas. [0134] Fig.11 illustrates the success rate in changing the network’s prediction when extremal deformations are applied on individual substructures of the brain. In Fig.11, the success rate in producing such a deformation that changes the DNN’s prediction from NC to AD is presented. Surprisingly, few of these deformations manage to change the prediction initially made by the neural network. This experiment suggests that the decision taken by a lightweight CNN stems from a combination of different image features that are spatially disjoint. The disclosed method highlights an asymmetry in the decision made by the DNN which might reflect an asymmetry in atrophy due to AD or a bias in the model. Conclusion [0135] This disclosure presents a post-hoc interpretability method based on the construction of smooth deformations of the input image data. To obtain localised and anatomically meaningful saliency maps, a counterfactual vector field norm is used. Going further, an alternative visualisation technique is proposed that allows monitoring the voxel changes in intensity between the original and counterfactual image. Finally, stability through local image deformation of the revered CNN architecture is established by a slight modification to the disclosed method. [0136] It will be appreciated by persons skilled in the art that numerous variations and/or modifications may be made to the above-described embodiments, without departing from the broad general scope of the present disclosure. The present embodiments are, therefore, to be considered in all respects as illustrative and not restrictive.