Upscaling the resolution of an image

Abstract

A computer-implemented method for upscaling the resolution of an image which uses a co-registered pair of images comprises a first image having a first resolution and a second image having a second resolution that is lower than the first resolution. The method determines local or global relationships between the images, and generates an upscaled version of the second image (a super-resolved image) using the local or global relationships and the pixel data of the first image. Embodiments of the invention may be applied to any imaging modalities, for example the first image and the second image may be derived from two different imaging modalities (e.g., images based on different wavelengths). For example, the first image may comprise a high-resolution intensity image and the second image may comprise a lower-resolution fluorescence lifetime imaging microscopy image.

Description

[0001] Upscaling the Resolution of an Image

[0002] Field of the Invention

[0003] The present invention generally relates to methods and apparatus for upscaling a low-resolution image based on correlations between the low-resolution image and a high-resolution image of a different modality.

[0004] Background

[0005] Computational super-resolution (SR) provides a route to obtain high image resolution from few physical samples. SR relies on solving an under-determined system, wherein there are M measurements, and the target data / signal needs to be reconstructed with N > M data points. Multiple flavours of SR exist, which are generally either interpolation, reconstruction- (inverse retrieval) or example- (learning) based. Interpolation is the simplest form of upsampling, encompassing a large family of methods for connecting datapoints with some curve [1], Although these methods increase the number of pixels in an image, they do not add new information. Reconstruction-based modelling instead manipulates the detection scheme to optically encode information about the high-resolution target pixels in fewer measurements. Manipulating the detection scheme typically requires hardware modifications, which can be cumbersome since they do not use conventional imaging techniques. As such, these modifications may be more difficult to implement and use with existing setups. Lastly, example-based schemes rely on computation and pattern recognition to upsample images in a data driven manner. However, learning-based schemes rely heavily on large, diverse training datasets. For natural images, there is a plethora of such data, but for niches in medical imaging, this can pose a bottleneck. Furthermore, example-based schemes are prone to out-of-distribution hallucination - estimating high-resolution features from the training set that are vastly different from the true sample features, yet which match the low-resolution measurements when downsampled.

[0006] An example situation where such super-resolution methods may be applied is fluorescence lifetime imaging microscopy (FLIM), which finds extensive applications in biological studies where the lifetimes of fluorophores can be used as indicators of the cell micro-environment. The generation of fluorescence lifetime images with time-resolved methods typically requires recovery of a full temporal decay curve for a suitable number of spatial points of the sample. Lifetime is then recovered with, for example, least squares deconvolution [2], Time correlated single photon counting (TCSPC) has emerged as one of the more common protocols for FLIM, in particular applied to laser scanning confocal microscopy. However, the scanning method presents drawbacks, such as the lack of instantaneous complete field-of-view information, or long acquisition times which are incompatible with the rapid intracellular dynamics of living cells.

[0007] The advancement of Single Photon Avalanche Diode (SPAD) array technology that incorporate TCSPC at each pixel, has paved way for wide-field FLIM. However, SPAD arrays typically have low-pixel counts and / or low fill-factor that impacts spatial resolution. In any case, increasing the native resolution TCSPC data also runs into data transfer / storage issues related to the 3-dimensional (2 spatial dimensions and 1 temporal dimension) nature of the data. Super-resolving the output of the SPAD array is useful to help address these problems, but the super-resolution methods themselves have the drawbacks outlined above.

[0008] The present invention has been devised in light of the above considerations.

[0009] Summary of the Invention

[0010] This invention introduces new methods and apparatus for upscaling the resolution of an image (which may be described herein as providing a super-resolved image). The goal of such an upscaling process is to produce an upscaled image which reproduces a ground-truth image, wherein the ground truth image represents a hypothetical image that would be taken by a high-resolution imaging device, as accurately as possible. At its most general, this is achieved by obtaining a co-registered pair of images, wherein the pair of images comprises a first image having a first resolution and a second image having a second resolution that is lower than the first resolution; determining local or global relationships between the images; and generating an upscaled version of the second image (a super-resolved image) using the local or global relationships and the pixel data of the first image. As will be apparent, the first and second images may be given any useful name, though herein the first image will be referred to as the reference image, and the second image will be referred to as a low-resolution image. Embodiments of the present invention may be applied to any imaging modalities, for example the first image and the second image may be derived from two different imaging modalities (e.g., images based on different wavelengths, and / or different data (e.g., obtained by a different procedure or method)). For example, the first image may comprise a high-resolution intensity image and the second image may comprise a lower-resolution fluorescence lifetime imaging microscopy (FLIM) image (wherein the ground-truth image is a FLIM image having the same resolution as the high-resolution intensity image). In another example, the first image may comprise a high-resolution intensity image at wavelengths in the visible spectrum and the second image may comprise a lower-resolution intensity image at infrared wavelengths (wherein the ground-truth image is an infrared image which may have the same resolution as the high-resolution intensity image).

[0011] The methods described herein may also be appropriate for up-scaling FLIM images, such as those produced from time-correlated single-photon capture (TCSPC) datasets obtained both from scanning confocal microscopy as well as single photon avalanche photodiode (SPAD) arrays. SPAD pixels generally have low-pixel counts and / or low fill-factor. The reference image may be produced by imaging equipment with a greater maximal resolution, which may be cheaper to produce, for example a CMOS sensor.

[0012] In a first aspect of the present invention, there is provided a computer-implemented method for upscaling the resolution of an image, wherein the method comprises: obtaining a co-registered pair of images, wherein the pair of images comprises a first image having a first resolution and a second image having a second resolution that is lower than the first resolution; determining a plurality of local correlation functions, wherein each local correlation function is associated with a respective sub-region of the second image and represents a correlation between pixel data of the respective sub-region with corresponding pixel data of the first image (in particular, spatially corresponding pixel data of the first image); and generating an upscaled version of the second image using the plurality of local correlation functions and the pixel data of the first image. The upscaled version of the second image many also include pixel-data of the second image. In this way, the present invention may produce greater cost efficiency, and / or higher resolution, images than other known imaging and / or upsamspling techniques. For example, the equipment used for producing the low-resolution image may have limited maximum resolution, or increasing the resolution may be prohibitively costly. The first image and the second image may be taken using different sensing modalities, for example representing different data types (e.g., intensity data of different wavelengths). By being co- registered, it may be understood that the first and second images represent the same, or an overlapping, field of view, such that pixels in the second image may spatially correspond to pixels in the first image. In examples, the upscaled version of the second image may have the same resolution as the first image.

[0013] Optionally, in a method according to a first aspect, generating an upscaled version of the second image may further comprise producing predicted pixel-data by applying each of the plurality of correlation functions to respective sub-regions of the first image, wherein generating an upscaled version of the second image is based on the predicted pixel-data. Generating an upscaled version of the second image may comprise using the predicted pixel-data to directly construct the upscaled version of the second image, or the predicted pixel-data may be used to construct the intermediate image, for example. Applying each of the plurality of correlation functions to respective sub-regions of the first image may comprise selecting a respective sub-region of the first image and producing predicted pixel-data for each pixel of the first image without a corresponding pixel in the second image, by applying a correlation function corresponding with the respective sub-region.

[0014] Optionally, a method according to a first aspect may further comprise determining a predicted error associated with the upscaling method. At least one local correlation function may be determined, the local correlation function being associated with a respective sub-region of the second image having a central pixel omitted. The local correlation function represents a correlation between pixel data of the respective sub-region with corresponding pixel data of the first image (in particular, spatially corresponding pixel data of the first image). The local correlation function may be used to predict pixel-data corresponding to the omitted central pixel of the sub-region. The method may include a step of comparing the predicted pixeldata with the corresponding sampled pixel-data of the second image to provide a statistical upper bound on the error which may be provided by the method of the first aspect of the present invention.

[0015] Optionally, in a method according to a first aspect, each respective sub-region may have a size of between 2 x 2 - 10 x 10 pixels of the second image. For example, the sub-regions may have a size of 5 x 5 pixels in the second image since this has been found empirically to provide good results over a range of upsampling factors.

[0016] According to a second aspect of the present invention, there is provided a method for upscaling the resolution of an image, the method comprising: obtaining a co-registered pair of images, wherein the pair of images comprises a first image having a first resolution and a second image having a second resolution that is lower than the first resolution; producing training data by obtaining a first plurality of patches in the first image, and, for each patch, obtaining a pixel in the second image corresponding to the centre of the patch; training an algorithm using the training data, the trained algorithm being configured to reproduce or predict pixel-data of the pixel in the second image corresponding to the centre of any patch of the first image; and generating an upscaled version of the second image based on the trained algorithm and pixel-data from the first image. Optionally, generating an upscaled version of the second image may further comprise producing predicted pixel-data, wherein predicted pixel-data is produced by identifying a second plurality of patches of the first image, each patch being centred on a pixel of the first image without a corresponding pixel in the second image; and supplying the second plurality of patches to the trained algorithm to predict pixel-data of a pixel in the second image corresponding to the centre of each patch of the second plurality of patches; and wherein generating an upscaled version of the second image is based on the predicted pixel-data. The upscaled version of the second image may also include pixel-data of the second image. In this way, the present invention may produce greater cost efficiency, and / or higher resolution, images than other known imaging and / or upsamspling techniques. For example, the equipment used for producing the second image may have limited maximum resolution, or increasing the resolution may be prohibitively costly. The first image and the second image may be taken using different sensing modalities, for example representing different data types (e.g., intensity data of different wavelengths). By being co- registered, it may be understood that the first and second images represent the same, or an overlapping, field of view, such that pixels in the second image may spatially correspond to pixels in the first image. In examples, the upscaled version of the second image may have the same resolution as the first image.

[0017] According to a third aspect of the present invention, there is provided a method for upscaling the resolution of an image, wherein the method combines steps of the methods described above with respect to the first aspect and the second aspect of the present invention. Specifically, the method according to the third aspect of the present invention comprises: obtaining a co-registered pair of images of different modalities, wherein the pair of images comprises a first image having a first resolution and a second image having a second resolution that is lower than the first resolution; determining a plurality of local correlation functions, wherein each local correlation function is associated with a respective sub-region of the second image and represents a correlation between pixel data of the respective sub-region with corresponding pixel data of the first image; producing training data by: obtaining a first plurality of patches in the first image; and, for each patch: obtaining a pixel in the second image corresponding to the centre of the patch; training an algorithm using the training data, the trained algorithm being configured to reproduce or predict pixel-data of the pixel in the second image corresponding to the centre of any patch of the first image; and generating an upscaled version of the second image using: the plurality of local correlation functions; the trained algorithm; and the pixel data of the first image.

[0018] Optionally, the algorithm may be a deep neural network. Such algorithms are well adapted for image processing, and can be either supervised, semi-supervised or unsupervised. Beneficially, this means that a method of the invention may be performed entirely without user input. By using a deep neural network, the method combines strengths of convolutional and fully connected neural layers. For example, convolutional layers at the start of the architecture may be good at hierarchical feature extraction (extracting shapes in a translation invariant manner), whereas fully connected layers can learn arbitrary non-linear algorithms to map these features onto the low-resolution image.

[0019] Optionally, producing the training data may further comprise rotating and / or reflecting the plurality of patches. Augmenting the training data produces a larger training data set, which leads to a more accurate trained algorithm. As such, the upscaling accuracy of the method is improved.

[0020] Optionally, a size of the patches is pre-determined based on features in the first image. This may be done either by a user or as a step in the computer-implemented method. Firstly, a morphological feature may be identified - this may be identified from the first or second images, or from knowledge of the image subject matter. The patches may then be determined from a size of the identified feature. Ideally a patch is chosen large enough to encompass this feature. In this manner, the algorithm may recognise features of the first image, leading to more accurate results.

[0021] Optionally, a method may further comprise determining a predicted error associated with the upscaling method. The method may include omitting a portion of the training data set used to train the algorithm. The omitted data may be used as a validation set by predicting pixel data of the validation set. The predicted pixel data produced by the trained algorithm may then be compared with the corresponding pixel-data of the low-resolution image. This provides a statistical upper bound on the upscaling error of a method according to a second aspect of the invention. Optionally, the validation set may be re-introduced into the training data and the algorithm may be retrained. This maintains the accuracy of the algorithm while providing a useful upper bound on the error of the method.

[0022] Optionally, a statistical upper bound on the error of a method may be used in further calculation. The upper bound may be compared to a threshold error. If the upper bound is greater than a threshold error, the method may include determining that an error level is too high, and alerting a user in response. The method may therefore alert the user when two provided images are not suitable for upscaling. For example, this may occur if there are not sufficient global or local correlations between either image.

[0023] Optionally, generating an upscaled version of the second image may comprise: constructing an intermediate image; optimizing the intermediate image based on the intermediate image, the first image, and at least one of the local correlation functions or the trained algorithm by varying the intermediate image; and outputting an optimized intermediate image as the upscaled version of the second image.

[0024] Optionally, the intermediate image may be optimised based also on the second image. In this way the method may further improve the accuracy of image upscaling, because the optimised image may be more consistent with the sampled data in the second image (the second image being interpreted as a sparsely sampled version of the ground truth image).

[0025] Optionally, the intermediate image may be optimised based also on the total variation of the intermediate image. In this way, the method may produce images having a lower total variation or less excessive or spurious detail. Since the ground truth image may have a low total-variation, the ground truth image may be better reconstructed by such a method. Optionally, optimizing the intermediate image may comprise constructing a cost function based on the intermediate image, the first image and the local correlation functions; and minimizing the cost function by varying the intermediate image.

[0026] Optionally, one or more terms in the cost function may be weighted, for example to achieve optimal results. In particular, at least one of the first image or second image may be weighted. That is, the optimisation may produce images which are closer to the priors, or closer to the sampled data by adjusting the weightings of the first image and / or the second image accordingly. In this way, the method may further improve on the upscaled image which is generated based on the priors. An optimal weighting factor may vary for the types of images used with the method, Therefore the weighting factor may be optimised to yield best results through empirical testing.

[0027] Optionally, obtaining a co-registered pair of images may comprise obtaining a pair of un-registered images; and co-registering the pair of un-registered images. In this manner the method may be adapted for use with images that have not been co-registered. This may be the case when the images have not been taken by the same sensor or by the same detector. Co-registering the images may include applying a transformation to the images, such as rotating or scaling the image, so that corresponding points in either image correspond to the same points of the FOV.

[0028] Optionally, the method may be used for upscaling images using an intensity image as a reference image. Typically, intensity images are cheaper to produce and may have greater maximal resolutions than other imaging techniques. Therefore, upscaling images using intensity images as a reference may be cheaper or produce higher resolutions than other known techniques. Likewise, the low-resolution image may be a lifetime image. It is known that lifetime images require expensive imaging equipment to produce and may not have as high resolutions as other images. Therefore, these images may be particularly suitable for the upscaling methods described herein.

[0029] Optionally, obtaining a co-registered pair of images may comprise: obtaining data that is not in the form of a single image; and computing, from the obtained data, the second image. This makes the method of image upscaling applicable to wider data types - especially those that are not in the form of images but may be represented as such. In one example, a lifetime image may be produced from a FLIM datacube.

[0030] Optionally, the data may obtained using a single-photon avalanche diode array. A single-photon avalanche diode array generally does not produce particularly high resolution images, and is therefore a good candidate for upscaling. However, this kind of data is not always produced in the form of an image, and therefore may need to be turned into an image before it can be upscaled.

[0031] In a fourth aspect of the invention, there is provided an imaging apparatus for capturing co-registered images and upscaling the resolution of an image comprising an objective lens configured to direct light from a subject; a beam-splitter configured to receive light from the objective lens and direct it to: reference imaging equipment configured to generate a reference image of the subject; low-resolution imaging equipment configured to generate a low-resolution image of the subject; and a processor configured to execute a computer-implemented method according to a first or second aspect. Such an apparatus may produce high-resolution images more cheaply or accurately than other known imaging techniques.

[0032] Likewise, by using a beamsplitter, the apparatus may be configured to generate a first and second image that are already co- registered.

[0033] In a fifth aspect of the invention, there is provided a computer system comprising a processor configured to perform a method according to a first or second aspect of the invention.

[0034] In a sixth aspect of the invention, there is provided a non-transitive computer-readable storage medium containing instructions which, when executed by a computer, cause the computer to carry out a method according to a first or second aspect of the invention.

[0035] The invention includes the combination of the aspects and preferred features described except where such a combination is clearly impermissible or expressly avoided.

[0036] Summary of the Figures

[0037] Embodiments and experiments illustrating the principles of the invention will now be discussed with reference to the accompanying figures in which:

[0038] Figure 1 is a flow-chart of a computer-implemented method according to an embodiment of the present invention;

[0039] Figure 2 shows a series of images representing an example application of the method shown in Figure 1 ;

[0040] Figure 3 is a flow-chart of a computer-implemented method according to a second embodiment of the present invention;

[0041] Figure 4 shows a series of images representing an example application of the method shown in Figure 3;

[0042] Figure 5 is an example of an upscaled version of a low-resolution image generated according to an embodiment of the present invention, and a comparison with prior-art upscaling methods;

[0043] Figure 6 is further example of an upscaled version of a low-resolution image generated according to an embodiment of the present invention, and a comparison with prior-art upscaling methods;

[0044] Figure 7 is a schematic diagram showing an apparatus suitable for enacting an embodiment of the present invention;

[0045] Figure 8 is an example of an upscaled version of a low-resolution thermal image generated according to an embodiment of the present invention, and a comparison with a prior-art upscaling method; and

[0046] Figure 9 depicts an example computing system. Detailed Description of the Invention

[0047] Aspects and embodiments of the present invention will now be discussed with reference to the accompanying figures. Further aspects and embodiments will be apparent to those skilled in the art. All documents mentioned in this text are incorporated herein by reference.

[0048] Disclosed herein is a data-fusion method for implementing computational super-resolution (SR), termed ‘single sample image fusion upsampling’ (SiSIFUS). SiSIFUS leverages statistical correlations between two measurements or images to evaluate a data-driven prior that informs an inverse retrieval algorithm about the non-sampled points of the scene. This approach may be used, for example, to super-resolve fluorescence lifetime images, employing data from two separate sensing modalities: high spatial resolution intensity images and low spatial resolution but high temporal resolution time correlated single photon counting (TCSPC) images. However, it will be appreciated that such methods may also be applied to other datasets, such as pairs of images which display different data types (e.g., the two images may be derived from two different imaging modalities, wherein the pixel data of the modalities may be correlated with one another) at different resolutions. The disclosed methods provide ‘single sample’ priors, i.e. functions or neural networks that are fitted on data pairs from the same sample. These functions may then be used to generate a super-resolved image. Optionally, this super-resolved image may be optimised (e.g., by constructing and minimising a cost function based on the prior) to produce an even more accurate image. In this way, the super-resolved image which is obtained may be an accurate approximation of a ground-truth image (i.e., the hypothetical image which would have been obtained from a high-resolution imager). Unlike existing methods, the disclosed methods produce explicit functional mappings between the first and second images.

[0049] Two types of SiSIFUS are disclosed. According to a first aspect of the present invention, which may be referred to herein as ‘local SiSIFUS’, a super-resolved image is derived based on local correlations which may exist between pixels in a low-resolution image and corresponding pixels in a high-resolution image. According to a second aspect of the present invention, which may be referred to herein as ‘global SiSIFUS’, a super-resolved image may be derived from global morphological signatures in the image.

[0050] For example, one embodiment uses a neural network to predict fluorescence lifetime from intensity patches, finding global morphological similarities in the intensity image that may best explain the available fluorescence lifetime imaging microscopy (FLIM) samples. SiSIFUS is effective on TCSPC datasets obtained both from scanning confocal microscopy as well as SPAD arrays.

[0051] The disclosed methods may be benchmarked against bilinear interpolation approaches that can be achieved on a FLIM datacube without data fusion, i.e. interpolation in FLIM space (FS), and interpolation in photon space (PS).

[0052] Fig. 1 is a flow-chart showing a local SiSIFUS method 100 according to an embodiment of the present invention. In particular, the method 100 may be computer-implemented, for example by being provided as a set of processor-executable instructions. Fig. 2 shows a series of images representing an example application of the method shown in Fig. 1. At a first step, a co-registered pair of images is obtained 115. Specifically, the co-registered pair of images comprises a reference image and a low-resolution image. Obtaining a co-registered pair of images may be achieved by accessing a storage medium to retrieve the pair of images, or may comprise a separate step of capturing first and second images. The reference image and low-resolution image may be captured using different imaging modalities (e.g., the images may be captured with different equipment and / or with different techniques).

[0053] Obtaining a pair of co-registered images 115 may involve of obtaining a pair of un-registered images 110 and co-registering the images 120. Un-registered images are images which have not been co-registered. In some examples, images may be co-registered by virtue of how they are captured and may not require an additional co-registering step. Both images may be taken of the same subject, and share a field-of- view (FOV). The pixel-data of pixels in the first and second images are based on measurements taken at a point of this FOV. By co-registered it will be understood that for each pixel in either of the first image or the second image, a corresponding point in the other image can be defined which nominally refers to a measurement taken at the same point of the field-of-view. In many cases, a pixel will not have a corresponding pixel in the other image (e.g., due to the different resolutions of the first image and the second image). Pixels in either image having pixel-data based on measurements taken at the same point of the FOV shall be referred to as registered pixels. Co-registering the images 120 may involve applying a transformation to either image, such as rotating or scaling or cropping the image, so that corresponding points in either image correspond to the same points of the FOV.

[0054] In some embodiments, obtaining a co-registered pair of images may comprise a pre-processing step. This step involves obtaining data that is not in the form of a single image (e.g a FLIM datacube), and producing a low-resolution image (e.g. a lifetime image). In one embodiment, data of the form of a FLIM datacube is pre-processed into a lifetime image which may be used as the low-resolution image to be upscaled.

[0055] Many forms of pre-processing may be considered, depending on the type of image that is to be upscaled. Examples of such pre-processing include obtaining a depth image from a time-resolved sonogram with a maximum projection; obtaining a range-doppler-angle map from a frequency modulated continuous wave (FMCW) radar using Fourier transformations; or obtaining a 3D x-ray image from a computed tomograph (CT) scan with the Radon transform. Other processes involving obtaining data that is not in the form of a single image and producing a low-resolution image will be apparent to those skilled in the art. In one example, a fluorescence lifetime imaging microscopy (FLIM) datacube may be converted into a lifetime image. This lifetime image may then be upscaled directly.

[0056] The resolution of the low-resolution image is lower than that of the reference image. The difference in resolution may be defined by an integer, for example the reference image may have a resolution 2x greater than that of the low-resolution integer. As a result, each of the pixels of the low-resolution image may have a corresponding pixel in the reference image.

[0057] By way of example, Figs. 2A depicts an overlay 410 wherein a low-resolution image is overlaid onto a high-resolution reference image. Each point of the low-resolution image has a corresponding point in the reference image, although not every pixel in the reference image has a corresponding pixel in the low- resolution image. Registering the points of the high-resolution image to the points of the low-resolution image may involve registering each point of the high-resolution image with a nearest pixel of the low- resolution image. In this way, the overlay 410 shows a pair of co- registered images.

[0058] The local SiSIFUS method 100 exploits local relationships that may exist between pixel-data of the low- resolution image and the corresponding pixel-data of the reference image. After obtaining a co-registered pair of images, the method 100 includes a step of determining a plurality of local correlation functions 130, wherein each local correlation function is associated with a respective sub-region of the low- resolution image and represents a correlation between pixel data of the respective sub-region of the low- resolution image with corresponding pixel data of the reference image.

[0059] Determining the local correlation functions 130 may include a step of identifying a plurality of sub-regions in the low-resolution image. A sub-region is a subset of neighbouring pixels within an image and may divide the low-resolution image into distinct sub-regions, or the sub-regions may overlap. For example, each sub-region of the plurality of subregion may be a square of 2 x 2 - 10 x 10 pixels in the low- resolution image. For each sub-region in the low-resolution image a corresponding sub-region in the reference image is defined (for example, a sub-region defining the same field of view), such that every registered pixel in the sub-region of the low-resolution image may have a corresponding pixel in the subregion of the reference image.

[0060] The dimensions of the sub-regions may be pre-determined. The dimensions should be chosen so as to ensure that a sub-region captures morphological features of interest when upscaling. The method may include determining a dimension of sub-region to be used, by first identifying a size of a morphological feature of the image. Alternatively, the method may include generating a plurality of upscaled images using sub-regions of different dimensions for each image. The upscaled version of the low-resolution image may be output as the upscaled image having the lowest error (for methods which involve calculating a statistical upper bound on the error of the method). This may not be used in all embodiments, since the mean average error of upscaling does not necessarily reflect human perception as well as other metrics such as learned perceptual image patch similarity (LPIPS).

[0061] Figure 2B depicts a subregion (or window) 420 of 5x5 pixels in the low-resolution image overlaid onto the corresponding sub-region of the reference image, which is 10x10 pixels. In this example, the method is thereby configured to provide 2x upscaling of the low-resolution image. From this, it will be appreciated that all pixels of the low-resolution image have a corresponding pixel in the reference image, but so do the areas between pixels of the low-resolution image. The method 100 may therefore use local correlation functions to fill the areas between pixels of the low-resolution image, in order to provide a super-resolved image directly, or to provide an intermediate image which can be used as the basis for optimisation which is described in more detail below.

[0062] Each correlation function is associated with a sub-region of the low-resolution image, Each and represents the relationship between pixel-data of a pixel of the low-resolution image and pixel-data of its corresponding registered pixel in the reference image, over the associated subregion. There may be many registered pixels, so determining each correlation function may involve determining a best fit for the pixel-data of registered pixels in this subregion. Each correlation function may be determined such that it may act on any input pixel-data, but it is applied in particular to pixels of the reference image which do not directly correspond with pixels of the low-resolution image. Applying the correlation function to pixel-data of the reference image results in new pixel-data, which may be used as a prediction or estimation of pixeldata of the ground-truth image (that is, a prediction or estimate of pixel-data associated with the areas between pixels of the low-resolution image).

[0063] Figure 2B depicts a graph 430 showing data from the sub-region 420 shown above. On the y-axis, the graph 430 plots pixel-data of the low-resolution image. On the x-axis the graph 430 plots the corresponding pixel-data of the reference image. Shown on the graph 430 are three possible correlation functions to ‘best-fit’ the data, which in this example include a linear interpolation, a cubic spline, and a radial basis function (RBF) gaussian process (GP). However, it will be appreciated that any correlation function which fits the data will be appropriate.

[0064] The correlation functions may be any suitable transformation, for instance: Gaussian Processes with RBF kernels, B-spline fitting, or Interpolation (Nearest, Linear, Cubic). Mean upsampling metrics over a validation set for multiple upsampling factors (e.g. 2x, 4x, 8x and 16x) may be used to test the performance of each transformation. In some examples, a metric may be used, for example as produced by leave-one-out cross-validation (LOOCV), to estimate the quality of upsampling which results from a given method.

[0065] The correlation functions may be used in a step of predicting, from the reference image, pixel-data of the ground-truth image in order to construct an intermediate image 140 or, in some examples, directly output this prediction as an upscaled version of the second image. Each sub-region of the reference image has a corresponding sub-region of the low-resolution image, and is thereby associated with a respective correlation function. Therefore, each respective correlation function may be applied to a respective subregion of the reference image in order to output (or produce) predicted pixel-data. This may occur across the whole respective sub-region, or just a subset of pixels within this sub-region. Specifically, each respective correlation function is applied locally on a respective sub-region of the reference image, which in turn corresponds to the sub-region of the low-resolution image associated with the correlation function.

[0066] The correlation functions then output (or produce) predicted pixel-data. The correlation functions may predict pixel-data of the ground-truth image with more accuracy than prior art methods (e.g., interpolation of the low-resolution image). As the correlation functions may better fit the correlation between pixel-data of the reference and low-resolution images, an upscaled image produced using these correlation functions may have improved accuracy compared to prior art methods.

[0067] Determining the plurality of correlation functions may include omitting a central pixel of a sub-region in the low-resolution image (this process may be known as leave-one-out cross-validation, LOOCV). A correlation may be determined between pixel data of remaining pixels in the sub-region of the low- resolution image and the corresponding pixels of the reference image. The method may include predicting, from the reference image and correlation functions, pixel-data of the ground-truth image, as described above. Comparing the pixel-data predicted with and without omitting a central pixel (e.g., by calculating mean-absolute-error) may provide a statistical upper bound on the error expected from the method 100 since the windows contain pixel-data for one fewer pixel and the central pixel is maximally distant from the other pixels in the window. Alternatively, the correlation functions produced omitting a central pixel may be used to generate the upscaled image without calculating an error.

[0068] Optionally, this upper bound on the error may be used in further calculation. The upper bound may be compared to a threshold error. If the upper bound is greater than a threshold error, the method may include determining that an error level is too high and alerting a user in response. A high threshold error may indicate that there is insufficient correlation between the two images to perform accurate upscaling. The alert provided to the user may include an indication that the two images obtained are not suitable for the disclosed methods.

[0069] Fig. 2C depicts a second overlay 440 wherein a low resolution image is overlaid onto a high resolution reference image, similar to the overlay 410. However, in the second overlay 440, the central pixel of the low resolution image is removed in order to facilitate an error calculation by LOOCV as described above.

[0070] The local SiSIFUS method 100 may, in certain embodiments, directly generate an upscaled image based on the local correlations. For example, pixel data may be predicted and directly used to generate an upscaled image. Generating an upscaled image based on a reference image introduces more information than methods based on a low-resolution image alone. However, in other embodiments such as that shown in Fig. 1 , the local correlations, with pixel-data of the reference image, are used to generate a ‘prior’ or intermediate image 140 which is then used as the input of an optimisation process to generate an improved upscaled image. The method 100 thereby produces upscaled images with greater accuracy to a ground truth image.

[0071] Generating an upscaled version of the low-resolution image 145 therefore comprises optimising an intermediate image 155, which has been generated based on the reference image and the local correlation functions. This super-resolved intermediate image may be output as an upscaled version of the low-resolution image, or it may be optimised for greater accuracy to a ground truth image.

[0072] Generating an upscaled version of the low-resolution image may be achieved iteratively, by modifying the intermediate image until eventually the solution converges on the upscaled image. An optimal upscaled image based on the reference image and either the trained algorithm or local correlation functions may not exist or be easily calculated. By converging on a solution rather than calculating directly, the method may produce an accurate upscaled image in a reasonable time. Empirically, an image iteratively optimised in this manner has produced more accurate or desirable images than the intermediate image without optimisation.

[0073] The intermediate image may also be optimised based on the low-resolution image, ensuring that the upscaled version of the low-resolution image does not deviate too significantly from the actual measured results which may improve the accuracy of the upscaling. Optimising the intermediate image may be achieved by constructing a cost function which is then minimised. This cost function may include terms upon which the optimisation should be based. By constructing a cost function, it makes it possible to introduce terms affecting the output of the optimisation step, and also to weight these terms to effect optimisation.

[0074] As shown in the embodiment of Fig. 1 , the method 100 optimises the intermediate image by generating a cost function 150. The value to be optimised may be a cost function based on the reference image, an intermediate image and the correlation functions, as described below. The cost function may then be minimised as a function of the intermediate image 160. After minimisation, the intermediate image may be taken to as the optimal approximation for the high-resolution image. Finally, the method 100 outputs the intermediate image as the upscaled version of the low-resolution image 170.

[0075] As detailed below in the mathematical treatment, the cost function includes a term based on the local correlation functions. This term may be the pixel-data predictions produced by the reference image and the local correlation functions. The cost function may also include a term based on pixel-data of the low- resolution image. This ensures that the upscaled version of the low-resolution image does not deviate too significantly from the actual measured results, which may improve the accuracy of the upscaling.

[0076] The cost function may include a term containing the 2D total variation (TV) of the intermediate image [3], This term may be weighted by a second weighting factor optimised to yield best results. The weighting factor may be similar for all low-resolution and reference images, or the weighting factor may be decided each time the method is applied. The 2D TV may be calculated using the intermediate image, and an operator representing the finite differences approximation of the horizontal and vertical image gradients. In most images there are smooth features that do not have strong variations in intensity, therefore having a low TV. On the other hand, images with excessive and possibly spurious detail may have a high TV. To avoid such excessive detail, the method may include the TV in the cost function. Constraining the cost function in this manner may result in the generated upscaled image having a lower TV and therefore have more desirable qualities than those generated without a dependence on TV.

[0077] The cost function may include a term based on the low-resolution image. The low-resolution image can be thought of as a sparsely sampled version of the ground-truth image which the method aims to reconstruct. This hypothetical sparse sampling may be represented as a matrix multiplication of the hypothetical ground-truth image with a matrix representing the sparse sampling operation, and pixel-data of the sparsely sampled low-resolution image may be a term of the cost function. Terms in the cost function may be weighted by a respective weighting factor. These weighting factors may be optimised to yield best results through empirical testing. The weighting factors may be similar for all low-resolution and reference images, or the weighting factors may be determined each time the method is applied. Choosing ideal weightings may comprise obtaining a set of test of images and optimising the weights for these test images by comparing the upscaled images generated with each set of weightings with a ground truth image. The weightings may then be fixed and used without further changes for new images. Methods involving weighting terms of a cost function allow certain terms (e.g., the reference or low- resolution images) to affect the optimised image. Where the weighting factors are chosen appropriately, this may lead to the generation of more accurate or desirable images.

[0078] Optimisation may be implemented using any commonly known algorithm for minimising a cost function. For example, this may be achieved by an alternating direction method of multipliers (ADMM) algorithm. Any suitable number of iterations may be used for the ADDM algorithm, ideally between 20 - 100 iterations. Optionally, the number of iterations may be chosen such that additional iterations produce insignificant change in the output, for example once the changes between subsequent iterations falls below a threshold level the ADMM algorithm may be stopped. The threshold level may be chosen to be any suitable value, for example between 0.01 and 1.00%, such as between 0.01 and 0.05%, depending on the degree of accuracy which is desired. In one example, a total of 20 iteration steps may be used for minimizing the cost function.

[0079] Fig. 5 shows a series of images (a)-(f) as an illustration of a local SiSIFUS method, as applied to a SKOV3 Rac1 Raichu validation sample. The fluorescence lifetime datacube may be recorded using a custom microscope setup, which uses a FLIMera 192xi28x326-pixel SPAD array (e.g., from HORIBA Scientific [4]). Simultaneously, a high spatial resolution sCMOS camera registers a 2048x2048-pixel image of the sample. The sCMOS image is spatially co-registered to match the FOV of the SPAD camera image. Fig. 5(a) shows a ground truth image 610, which in this example is a real high-resolution image by way of example only, which is down-sampled to provide a low-resolution image 620. In this way, the efficacy of the local SiSIFUS method may be demonstrated, but it should be appreciated that the present invention does not require a real ground-truth image to be provided in this manner.

[0080] A least squares deconvolution may be applied to fluorescence lifetime datacube at each pixel to obtain a ground truth lifetime image 610, TGT. This image is then sparsely sampled to provide a low-resolution lifetime image 620 TLR fortesting various upsampling factors. This image may be upsampled using a local SiSIFUS (e.g., the method 100 described above with respect to Fig. 1), using a 5 x 5 pixel window size. Comparing the upsampled image to the original TGT explicitly evaluates the quality of the reconstruction. An output of the local SiSIFUS method for 2x upsampling is shown in Fig. 5(f).

[0081] The quality of the method may be established numerically with learned perceptual image patch similarity (LPIPS), a metric that computes similarity between the feature spaces of two images passed through a deep neural network [5], A lower LPIPS score indicates better perceptual similarity. LPIPS is preferred over other metrics based e.g. on L1 or L2 norms as these are known to for example to be relatively insensitive to image blur, an important factor when considering microscopy images.

[0082] Local SiSIFUS for upsampling factors of 2x, 4x, 8x, and 16x may be compared against the LPIPS of bilinear interpolation in FLIM space (FS) and photon space (PS) as shown by the graph 660 in Fig. 5(c), which represents the LPIPS score by comparing the ground truth image 610 with the output of each method. In FS interpolation, the low-resolution fluorescence 3D datacube is first fitted at each pixel to provide a low-resolution 2D lifetime image, and then the low-resolution lifetime image is spatially interpolated to produce a final FS image 630 in Fig. 5(d). In PS interpolation, the low-resolution datacube is first spatially interpolated at each time-slice, and then fitted to the interpolated datacube to obtain a high-resolution lifetime image as a final PS image 640 in Fig. 5(e). This second approach is found to outperform FS interpolation in all instances. A local SiSIFUS method, however, successfully outperforms both forms of bilinear interpolation for all the upsampling factors, as shown by the graph 660. Experimentally, 2x local SiSIFUS may outperform 2x bilinear FS by 45%, and 2x bilinear PS by 26%.

[0083] Upsampling factors of 2x, 4x, 8x, and 16x are shown in graph 660. A first line 661 shows the LPIPS scores for each upsampling factor using a local SiSIFUS method; a second line 662 shows the LPIPS scores for each upsampling factor using a bilinear photon space method; and a third line 663 shows the LIPIPs scores for each upsampling factor using a bilinear FLIM space method.

[0084] The local correlation priors produced by the local SiSIFUS method rely on locally sampling a given intensity feature in lifetime space. In contrast, global SiSIFUS uses global morphological priors (GMPs), which exploit correlations between morphological features of the reference image and pixel-data of the low-resolution image. As a result, global SiSIFUS may produce more accurately upscaled images than methods which do not exploit morphological features of the reference image.

[0085] Fig. 3 is a flow-chart showing a global SiSIFUS method 200 according to another embodiment of the present invention. In particular, the method 200 may be computer-implemented, for example by being provided as a set of processor-executable instructions. Fig. 4 shows a series of images as an example of the global SiSIFUS method 200.

[0086] At a first step of the global SiSIFUS method 200, a co-registered pair of images is obtained 115, generally as described above with respect to the local SiSIFUS method 100 shown in Fig. 1 .

[0087] After obtaining a co-registered pair of images 115, the method 200 includes a step of producing training data 225. Each individual piece of training data comprises a patch in the reference image and a pixel in the low-resolution image corresponding to the centre of the patch. The training data therefore comprises a plurality of pieces of training data using a plurality of patches of the reference image, each with a corresponding pixel of the low-resolution image. An algorithm is trained on this training data as part of a training step 230.

[0088] By way of example, Fig. 4A depicts an overlay 510 wherein a low-resolution image is overlaid onto a reference image, demonstrating the correspondence between the pixels in either image. Each point of the low-resolution image has a corresponding point in the reference image, although not every pixel in the reference image has a corresponding pixel in the low-resolution image. A first plurality of patches 520 show patches of the reference image having central pixel corresponding to a pixel in the low-resolution image.

[0089] To deal with the relatively small amount of training data contained in a single sample pair of images, the method 200 may include augmenting the training set. For example, producing training data 225 may include a further step of data augmentation. Any suitable data augmentation technique may be used, for example reflecting and rotating the reference image patches to form an augmented data. This increases the training data eightfold, as each patch of training data has four rotational versions, which each rotational version having a mirrored or reflected version.

[0090] Fig. 4B shows a training set 530 consisting of eight images, wherein the training set 530 is based on a single patch which has been rotated and flipped. It will be appreciated that the full set of training data could be based on every patch taken from the overlay 510.

[0091] The method 200 may produce further training data by identifying neighbouring patches of the reference image to those already identified. This may include identifying a neighbourhood of (or, an area around) a registered pixel in the reference image, wherein the registered pixel corresponds with a pixel of the low- resolution image. For each pixel of the reference image in this neighbourhood, the computer may identify a patch with this pixel at the centre. There may be no corresponding pixel-data in the low-resolution image for this pixel. Instead, training data is formed from this patch together with the pixel-data of the neighbouring registered pixel in the low-resolution image. That is, the training set may be augmented by estimating the data of the patches neighbouring the sampled patches with the same data as their sampled neighbour. This dataset augmentation method assumes that within the pixel-pitch of the intensity image, the variable measured by the imaging modality for the low-resolution image changes slowly. Using this augmentation technique, a neighbourhood of size m gives m2more windows to train on. Combined with the previous dataset augmentation, these neighbourhoods yield 72x and 200x more samples for 3x3 and 9x9 upsampling, respectively. This form of blind labelling can increase the diversity of the training input set. However, it should be noted that the true measurement data of neighbouring patches is unknown, so the estimated labels have inherent uncertainty. Consequently, this form of augmentation may be more effective when the training set size is very small, e.g. upsampling factors of 8 and above.

[0092] However, it should be noted that the pixel-data for the neighbouring pixels as obtained by this estimation is not necessarily the same as the ground-truth image which would be captured if an imaging device taking the second image was a higher resolution, as the true data for the neighbouring pixels is unknown. This may therefore introduce some inherent inaccuracy, which may be offset to some degree by using a larger training data set. Using neighbouring patches may result in a more accurate trained algorithm for high upsampling factors of 10 and above, since the deviation of a central pixel from its neighbour is likely to be smaller the higher the resolution of the ground-truth image. Using a more accurate trained algorithm to generate the upscaled image may increase the accuracy of the upscaled image.

[0093] Fig. 4B shows an example of a second training set 430 which is obtained by considering the neighbourhood around the same patch as for the training set 420. These nearest windows provide additional image patches to use fortraining data, and it can be seen that in the neighbourhood of the central pixel, the intensity (from the reference image data) is generally the same as for the central pixel, and so the assumption is that the lifetime (measured by the low-resolution FLIM image) also does not vary within the pixel-pitch of the reference image. It will be appreciated that the full set of training data could be based on every patch taken from the overlay 510.

[0094] Aligning the patches of the reference image with pixel-data of the low-resolution image at the centre of the patch may exploit correlations between reference image morphology and the low-resolution image. Real-life subjects often contain multiple examples of the same features in the reference image which correspond to similar features in the ground-truth image, and therefore exploiting correlations between reference image morphology and the low-resolution image may provide more accurate predictions for the ground-truth image.

[0095] The algorithm which is trained with the training data is an algorithm configured to predict pixel-data of a registered pixel in the low-resolution image which corresponds to the centre of any patch (which may be an arbitrarily chosen patch) of the reference image. For example, the algorithm may be a deep neural network (DNN) with random initialisation, and not a model trained on data from other images. The trained algorithm may be trained exclusively on the training data obtained from the reference image patches and pixel-data of the centre pixel of the corresponding low-resolution image patch obtained from a single sample (i.e., a single pair of images representing a single field of view). This improves the accuracy of the trained algorithm and may reduce artificial hallucination since the trained algorithm does not introduce biases from other datasets.

[0096] Figure 4C depicts an example step of training an algorithm based on training data 230. When provided a patch 560 of the reference image, the algorithm (in this case a deep neural network) is trained to reproduce pixel-data (e.g., in this example a lifetime estimate of 0.28) of the registered pixel in the low- resolution image corresponding to the centre of the patch. The algorithm may then take in patches and predict pixel-data for the central pixel. Repeating this process for multiple neighbouring patches predicts pixel-data for a portion of the ground-truth image. For example, a set of patches with unsampled centres 570 may be input to the DNN to give corresponding lifetime estimates 580 forthose patches. This may be repeated across the whole of the reference image.

[0097] The trained algorithm may be used in a step of predicting, from the reference image, pixel-data of an upscaled, or super-resolved, image. Predicting pixel-data may be achieved by first identifying a second plurality of patches of the reference image. Each of these patches may be centred on a pixel of the reference image which does not have a corresponding pixel in the low-resolution image. Each of these patches may be supplied to the trained algorithm which may predict pixel-data corresponding to the central pixel of each patch.

[0098] Near the edges of the reference image there are pixels which are not at the centre of any patch of the reference image. As such, the method does not involve predicting pixel-data for these pixels.

[0099] Alternatively, linear interpolation may be used to generate pixel-data for these pixels. Where this may not be possible, the methods described herein may involve only upscaling a central portion of the low- resolution image.

[0100] The method may include a step of predicting pixel-data of a patch of the reference image centred on a pixel of the reference image with a corresponding pixel in the low-resolution image. Comparing the predicted pixel-data with the value of the corresponding pixel-data in the low-resolution image may provide an estimate of the error produced by the method. This may be implemented by removing part of the augmented data set described above, and not using it as part of the training data used to train the algorithm. This removed data may be used as a validation set to test the accuracy of the trained algorithm. Using the trained algorithm to predict pixel-data for patches of the validation set and comparing this predicted pixel-data to the corresponding pixel-data of the low-resolution image provides a statistical upper bound for the accuracy of the method.

[0101] Optionally, the removed data may be reintroduced to the training set and the trained algorithm may be retrained on the whole training set. The retrained algorithm is expected to have a smaller error than that estimated by the validation set, since the overall training set used to train the algorithm is larger.

[0102] Optionally, this upper bound on the error may be used in further calculation. The upper bound may be compared to a threshold error. If the upper bound is greater than a threshold error, the method may include determining that an error level is too high, and alerting a user in response.

[0103] The method 200 may be used with new subjects which have previously unseen morphologies, reference image morphology-low-resolution image correlations, and low-resolution image ranges. Using an algorithm specialised on a single image ensures that the algorithm does not carry over biases from some external training set. Generating an upscaled image based on this trained algorithm may have increased accuracy compared to prior art methods.

[0104] The global SiSIFUS method is designed to generalise to new samples which have previously unseen morphologies, data correlations, and data ranges. Importantly, the DNN is just an architecture with random initialisation and not a model trained on data from other images.

[0105] For example, the DNN trained on data shown in Fig. 4B is trained exclusively on intensity windows and lifetime values obtained from a single sample (i.e. a single microscope field of view). To predict the lifetime in a pixel that the SPAD did not sample, the intensity patch may be found which is centred on this unsampled pixel, this patch passed through the trained DNN, as shown in Fig. 4C. Therefore, this scheme ensures that the model does not carry over biases from some external training set. For example, the DNN may be trained on intensity patches 13 x 13 pixels large, with ADAM [6], using a batch size of 100 over 150 epochs with MAE as the training loss. The size of the intensity patches are chosen to include the important features of the reference image. For example, as shown in Fig. 4A, there are repeated spots in the intensity image. The patches, shown in Fig. 4B, have been chosen to be as small as possible (so as to have as many windows as possible from the same image) but large enough to include each of these spots.

[0106] As with the local SiSIFUS method 100, the global SiSIFUS method 200 may, in certain embodiments, directly generate an upscaled image, based on the global correlations. For example, pixel data may be predicted and directly used to generate an upscaled image. Generating an upscaled image based on a reference image introduces more information than just upscaling a low-resolution image alone. However, in the depicted embodiment of Fig. 3, the global correlations, with pixel-data of the reference image, may be used to generate a ‘prior’ or intermediate image 240 which is then used as the input of an optimisation process to generate an improved upscaled image. The method 200 therefore produces upscaled images with greater accuracy to a ground truth image. Generating an upscaled version of the low-resolution image 245 therefore comprises optimising an intermediate image 255, which has been generated based on the reference image and the trained algorithm. The intermediate image may also be based on the low-resolution image, ensuring that the upscaled version of the low-resolution image does not deviate too significantly from the actual measured results which may improve the accuracy of the upscaling.

[0107] Generating an upscaled version of the low-resolution image may be achieved iteratively, by modifying the intermediate image until eventually the solution converges on the upscaled image. An optimal upscaled image based on the reference image and either the trained algorithm or local correlation functions may not exist or be easily calculated. By converging on a solution rather than calculating directly, the method may produce an accurate upscaled image in a reasonable time. Iteratively optimising an intermediate image, using the trained algorithm or local correlation functions as a constraint, may produce more accurate images than those produced by either the trained algorithm or local correlation functions alone.

[0108] The method therefore comprises generating a cost function 250. The value to be optimised may be a cost function based on the reference image, an intermediate image and the trained algorithm, as described below. The cost function may then be minimised as a function of the intermediate image step 260. After minimisation, the intermediate image may be taken to as the optimal approximation for the high-resolution image. Finally, the method 200 therefore outputs the intermediate image as the upscaled version of the low-resolution image 270.

[0109] As detailed below in the mathematical treatment, the cost function includes a term based on the trained algorithm. This term may be the pixel-data predictions produced by the reference image and the trained algorithm. The cost function may also include a term based on pixel-data of the low-resolution image. This ensures that the upscaled version of the low-resolution image does not deviate too significantly from the actual measured results, which may improve the accuracy of the upscaling.

[0110] The cost function may include a term containing the 2D total variation (TV) of the intermediate image. This term may be weighted by a second weighting factor optimised to yield best results. The weighting factor may be similar for all low-resolution and reference images, or the weighting factor may be decided each time the method is applied. The 2D TV may be calculated using the intermediate image, and an operator representing the finite differences approximation of the horizontal and vertical image gradients. In most images there are smooth features that do not have strong variations in intensity, therefore having a low TV. On the other hand, images with excessive and possibly spurious detail may have a high TV. To avoid such excessive detail, the method may include the TV in the cost function. Constraining the cost function in this manner may result in the generated upscaled image having a lower TV and therefore have more desirable qualities than those generated without a dependence on TV.

[0111] The cost function may include a term based on the low-resolution image. The low-resolution image can be thought of as a sparsely sampled version of the ground-truth image which the method aims to reconstruct. This hypothetical sparse sampling may be represented as a matrix multiplication of the hypothetical ground-truth image with a matrix representing the sparse sampling operation, and pixel-data of the sparsely sampled low-resolution image may be a term of the cost function. Terms in the cost function may be weighted by a respective weighting factor. These weighting factors may be optimised to yield best results through empirical testing. The weighting factors may be similar for all low-resolution and reference images, or the weighting factors may be determined each time the method is applied. Choosing ideal weightings may comprise obtaining a set of test of images and optimising the weights for these test images by comparing the upscaled images generated with each set of weightings with a ground truth image. The weightings may then be fixed and used without further changes for new images. Methods involving weighting the cost function may lead to the generation of more accurate or desirable images.

[0112] Optimisation may be implemented using any commonly known algorithm for minimising a cost function. For example, this may be achieved by an alternating direction method of multipliers (ADMM) algorithm. Any suitable number of iterations may be used for the ADDM algorithm, ideally between20 - 100. Optionally, the number of iterations may be chosen such that additional iterations produce insignificant change in the output, for example once the changes between subsequent iterations falls below a threshold level the ADMM algorithm may be stopped. The threshold level may be chosen to be any suitable value, for example between 0.01 and 1.00%, such as between 0.01 and 0.05%, depending on the degree of accuracy which is desired. In one example, a total of 20 iteration steps may be used for minimizing the cost function.

[0113] Fig. 6 shows a series of images (a)-(f) as an illustration of a global SiSIFUS method for upsampling, as applied to a Convallaria sample dyed with Acridine Orange, by way of comparison with other upsampling methods. The measurement may be performed using the previously mentioned SPAD array and sCMOS based microscope, yielding a datacube of size 192 x 256 x 326. A spatially co-registered sCMOS intensity image of dimensions 192 x 256 is also used. Fig. 6(a) shows a ground truth image 710, which in this example is a real high-resolution image by way of example only, which is down-sampled to provide a low-resolution image 720. In this way, the efficacy of the global SiSIFUS method may be demonstrated, but it should be appreciated that the present invention does not require a real ground-truth image to be provided in this manner.

[0114] The ground truth image 710 obtained from the SPAD array is shown in Fig. 6(a), and the corresponding 8x sparsely sampled low resolution lifetime image 720 is shown in Figure 6(b). The graph 760 shows the upsampling LPIPS scores indicating that global SiSIFUS may perform as well as (for low upsampling) or better (for higher upsampling) than bilinear interpolation (either in FLIM space (FS) or photon space (PS), as described above).

[0115] Figs. 6(d)-(f) show the results of 8 x 8 upsampling of the three compared methods. Fig. 6(d) shows a final FS image 730; Fig. 6I shows a final PS image 740; and Fig. 6(f) shows a final global SiSIFUS image 750. The 8 x 8 upsampling results in 730, 740 and 750 illustrate that global SiSIFUS may recognise that globules tend to have higher lifetimes than cell walls, hence global SiSIFUS may maintain contrast between these structures more consistently than both modes of bilinear interpolation. It may be noted, however, that global SiSIFUS may miss certain hotspots in the ground truth lifetime image 710 as might be expected, i.e. very small, sharp features are hard to reconstruct. It is also noted that since the global SiSIFUS method trains a model from scratch each time, different upsampling instances are of varying quality. To ensure the results are representative, the DNN architecture was retrained 3 times with different random initialisations for each upsampling factor, and the outputs with median LPIPS were used for the plotted results.

[0116] Even with weak local spatial correlation, there may be a strong global trend linking intensity patches and lifetime (for example, cell membranes have high lifetimes, and vesicles have low life-times). This may be exploited by the DNN prior to transform the low-resolution image into a high-quality reconstruction.

[0117] In another embodiment of the present invention (discussed in the mathematical treatment), the methods of both Figs. 1 and 3 may be performed together, for example in parallel, to determine a plurality of local correlation functions and obtain a trained algorithm. The upscaled image may then be generated based on both the local correlation functions and the trained algorithm along with pixel data of the first image. The method may include directly generating an upscaled image, based on the local correlation functions, global correlations and pixel data of the first image. For example, pixel data may be predicted by both local and global methods (as described above), which may then be combined, for example by averaging the two predictions. Generating an upscaled image based on a reference image introduces more information than just upscaling a low-resolution image alone. Generating an upscaled image based on both the trained algorithm and the local correlation functions allows an image to be generated based on local correlations and morphological features, which may lead to a more accurately upscaled image than images generated based on either the local correlations or morphological features alone.

[0118] Optionally, the method may include generating a cost function based on the reference image, an intermediate image, trained algorithm and local correlation functions. This cost function may be optimised in the same manner as described in relation to the local and global methods. As with the local and global methods, the intermediate image may be optimised based also on the second image, the total variation of the intermediate image, or the intermediate image may be optimised by minimising a cost function based on the first image and second image, and wherein at least one of the first image or the second image is weighted (as described above in relation to the local and global methods).

[0119] Figure 7 depicts an apparatus 900 for capturing a reference image and a low-resolution image which may be upscaled according to methods of the present invention. Specifically, the apparatus 900 may be used to capture an intensity image of a sample 910 as a reference image, and a FLIM image of the sample 910 as a low-resolution image. In this way, a method according to an embodiment of the present invention (e..g, a method as described above with respect to Fig. 1 or Fig. 3) may be used to obtain a superresolved FLIM image of the sample 910.

[0120] The apparatus 900 comprises an objective lens 940 and a beamsplitter 970 configured to direct light from the sample 910 to reference imaging equipment 920 and low-resolution imaging equipment 930. Light from the sample 910 is focussed onto the reference imaging equipment 920 by a reference lens 925, and focussed onto the low-resolution imaging equipment 930 by a low-resolution lens 935. The low-resolution imaging equipment 930 and reference imaging equipment 920 are configured to capture low-resolution images and a reference images, respectively. A processor in the apparatus then executes a method according to a first or second aspect of the present invention (e.g., a local SiSIFUS method 100 or a global SiSIFUS method 200) on these images in order to obtain a super-resolved image. The arrangement of beam-splitter, reference lens and low-resolution lens may be such that the images obtained by the apparatus do not require co-registering.

[0121] A laser 950 is configured to emit light through the objective lens 940 and onto the sample 910, via a laser focussing lens 955 and a half-silvered mirror 965, exciting flourophores within the sample 910. The sample 910 may then be imaged by the low-resolution imaging equipment 930 to form a FLIM image. A mirror 960 is used to re-direct light from the sample 910 through the beamsplitter 970 to help ensure that the apparatus 900 is compact.

[0122] The SiSIFUS method may also been applied for upscaling thermal images, for example. Fig. 8 shows a series of images as an illustration of a local SiSIFUS method as applied to thermal imaging. Thermal images in the dataset have been acquired through a Flir Vue Pro camera and the visible images through a Sony Xperia camera. In general, thermal images map the infra red energy content of a scene via a heat sensor attached to lens which then is integrated with standard image-capture technology.

[0123] Keeping in line with the general principle of SiSIFUS, a test low resolution thermal image may be upsampled by leveraging the information present in a different image of the same scene but acquired through different modality, which in this case may be a standard visual image.

[0124] Fig. 8 (b) shows a low resolution thermal image 1020, obtained by downsampling (via binning 8x 8 regions) a high resolution ground truth image 1010 (shown in Fig. 8(a)). A co-registered high resolution (256 x 336) intensity image 1030 is shown in Fig. 8(c), and local correlations were established between small patches (3 x 3) of the low resolution thermal image 1020 and the corresponding patch in the high resolution intensity image 1030 which may then be used to construct an intermediate image. The intermediate image may be fed to an inverse retrieval pipeline that aims to get the upsampled verison 1050 of the thermal image 1020. Fig. 8(d) also shows a bilinearly interpolated upsampled image 1040 for the same upsampling factor.

[0125] Mathematical Treatment

[0126] A general mathematical treatment demonstrating the principles underlying the present invention is included below.

[0127] The forward model that describes the sensing of the sample information may involve two measurements: the low spatial resolution image TLR (m,n) where m, n denote the pixels of the low-spatial resolution sensor. The second measurement is the high spatial resolution intensity measurement, l(M,N) where M, N denote the pixels of the high-spatial resolution sensor. TLR can be considered as the sparsely sampled version of a high resolution image (THR(M,N)) that the method aims to reconstruct, which may be referred to herein as the ground-truth image.

[0128] For example, in one embodiment of the present invention, the low spatial resolution image may be a low resolution FLIM image. This may be denoted by r(m, n, t), where m, n denote the pixels of the low resolution sensor, and t denotes the different time bins (e.g., the time of arrival of individual photons). In other embodiments, the first image and the second image may be derived from other imaging modalities, such as different wavelengths or the like.

[0129] The low-resolution image formation from its high-resolution counterpart can be expressed in the matrix operator notation as:

[0130] TLR = ATHR where A represents a sparse sampling operation. The two datasets, l(M,N) and TLR (m,n), are fed to a prior-generating pipeline (a local SiSIFUS method and / or a global SiSIFUS method as outlined above), that provide the best first-order approximations Tip (via the local prior generation method) and Tgp(via the global prior generation method) for THR. Local and global correlations are combined using both the LCP (local correlation prior) and GP (global prior) together in the loss function that is minimised in the inverse retrieval. The recovery of the optimal high resolution lifetime image THR* can be modelled as an optimization problem where the cost function that we minimize is:

[0131] THR* = arg min C(THR), where

[0132] C(THR) = 0.5|| ATHR - TLR II2i + 0.5y|| THR - Tip ||2i + 0.5p|| THR - Tgp|| i + all D THR HI, subject to THR S 0

[0133] The first term in C(THR) ensures the data fidelity between the low-resolution measured image and the downsampled optimal high resolution solution in each iteration. Prior constraints on the target high resolution lifetime image is enforced through the second and the third terms, weighed by the factors y and p which are empirically optimised to yield best results. The fourth term is the L1 norm of the 2D total variation (TV) evaluated on the high-resolution lifetime image and weighed by a. The anisotropic form of the TV is considered, and so the operator D represents the finite differences approximation of the horizontal and vertical image gradients. The optimization may be implemented using the alternating direction method of multipliers (ADMM) algorithm. A total of 20 iteration steps may be used for minimizing the cost function - for steps beyond that iteration number there may be an insignificant change in the solution.

[0134] For ADMM the minimisation problem can be re-formulated as

[0135] THR* = argmimuR 0.51| ATHR - TLR II i + 0.5y|| THR - Tip || i + 0.5p|| THR - Tgp|| i + a||z ||i, subject to D THR - z = 0 and THR > 0. The Augmented Lagrangian for this problem can be written as: LP(THR, Z, y) - 0.5|| ATHR - TL ll22 + 0.5y|| THR - Tip ||22 +0.5p|| THR - Tgp||22 + a||z ||i + y(D THR - Z)

[0136] +0.5p|| D THR - Z||22

[0137] Where y is a Lagrangian multiplier (or the dual variable) and p is a penalty parameter. The ADMM approach involves jointly minimizing the Lagrangian over all the primal variables followed by the updates over the dual variables. The primal updates for the variables THR and z are given by:

[0138] THR (k+1) <— argmin THR 0.5|| ATHR - TLR ll22+0.5y|| THR - Tip ||22 + 0.5p|| THR - Tgp||22 + y(D THR - ZI<) +0.5p|| D THR - zi< ll22

[0139] Zk+i <- argmiriz a||z ||i + y(D THR^ - Z) + 0.5p|| D THR - Z ||22

[0140] The dual update is given by: yk+1 <— yk + p (D THR - z)

[0141] For the primal minimisation updates, a standard optimization technique can be used based on the fast iterative soft thresholding algorithm (FISTA). Each iteration of the ADMM hence comprises of 90 iterations of FISTA for the THR variable update.

[0142] Local SiSIFUS as applied to FLIM images relies on the direct correlations between lifetime and intensity in real-world biological samples. Whilst there is no analytical formula linking fluorescent intensity with lifetime and literature often explicitly underlines the independence of these two parameters when measuring the fluorescent molecules in isolation, the inventors have realised that these variables are often locally correlated in real biological samples. Such local correlations may be seen by eye, both in FLIM videos (e.g. calcium imaging), and single frame measurements (e.g. SNARE imaging). Moreover, there are scenarios where fluorescence lifetime and intensity are explicitly linked. For instance, in Forster Resonance Energy Transfer (FRET), intensity and lifetime are directly proportional, though fluorescent dye concentration acts as a covariate, heavily influencing intensity but not lifetime.

[0143] A local correlation prior (LCP), or intermediate image, is generated by finding mappings between intensity and lifetime on small neighbourhoods (windows) of intensity-lifetime pixel pairs. A SPAD camera may sparsely sample a field of view: for a downsampling factor of A x A, lifetime is sampled at every A x Athpixel of the intensity image. The lifetimes corresponding to in-between pixels is unknown and must be statistically approximated - by upscaling the low-resolution FLIM image. Fig. 2 illustrate an example of the local SiSIFUS method as applied to FLIM images. An example pair of TLR and I are shown in Fig. 2A by the overlay 410 where the sparse, low-resolution lifetime samples are overlayed on the intensity image. The field of view (FOV) may be divided into windows, and spatially corresponding intensity-lifetime samples may be obtained. In each window, the intensity and lifetime pairs may be vectorised and fitted with a function, f, of the form Tip = f(l) (Fig.4(b)). Thus, our lifetime estimate Tip for pixel (Ai + x, Aj + y), is:

[0144] TiP(Ai + x, Aj + y) = f IJ (I (Ai + x, Aj + y)) with samples i e {0, 1 , ..., m -1} and j e {0, 1 , ..., n -1}, and x > 0, y < A. Importantly, the functions fij are fitted locally, not globally.

[0145] In lieu of an analytical formula linking lifetime and intensity, the form of these local functions may be found empirically from the data. A set of different function forms may be tried, including Gaussian Processes with RBF kernels, B-spline fitting, and Interpolation (Nearest, Linear, Cubic). The best form may be found by comparing mean upsampling metrics over a validation set of four biological samples (shown in Fig. 4B) for four upsampling factors each (2x, 4x, 8x and 16x). Experimentally, linear interpolation yielded the best results, hence it is used in all results shown in this disclosure.

[0146] The best function to map fluorescent intensity onto fluorescence lifetime with local SiSIFUS may be found by study. Tested functions may follow a set of common schemes, ranging from B-splines (linear, quadratic and cubic), through regular interpolation (nearest, linear and cubic); and kriging (radial basis function Gaussian processes, taking the mean of the fit) [7],

[0147] These functions may be applied to 4 samples and 4 upsampling factors (2x ,4x ,8x, and 16x). These methods may be evaluated based on mean-absolute-error and learned perceptual image patch similarity (LPIPS) between the reconstruction and ground truth, averaged over these 4 x 4 scenarios.

[0148] Local SiSIFUS can be modified for leave-one-out cross-validation (LOOCV), as in Fig. 4C. For this, the central intensity-lifetime measurement in the window is omitted from the fitting points, and f is fitted to the rest of the intensity-lifetime points in the window. With f, the central pixel’s lifetime T can be estimated from its intensity. By repeating this process for all windows, a prediction TLOOC is obtained.

[0149] The LOOCV mean-absolute-error (MAE), ||T- TLOOCV ||I, may provide a statistical upper bound on the reconstruction error one can expect with local SiSIFUS , | |THR— THR||I , since the windows contain one fewer intensity-lifetime sample and also the central pixel is maximally distant from the other samples in the window. Crucially, LOOCV is an a-priori metric, estimating the quality of upsampling without knowing the high-resolution ground truth.

[0150] Leave-on-out cross-validation (LOOCV) may act as an a priori uncertainty upper bound. The LOOCV scheme acts analogously to regular upsampling, except the central intensity-lifetime pair is removed from the observed data, and instead that pixel’s lifetime is predicted as though it were unknown. This has two effects. First, the number of intensity-lifetime samples in the window decreases by 1 . Second, the distance of the test pixel from its neighbours increases compared to other test pixels. The central pixel in a window is 1 pixel pitch (PP) away from its horizontal and vertical neighbours, and ^2x PP from its diagonal neighbours, whereas normally, each intensity pixel is at most 0.5 2 PP away from its nearest neighbour, and most are much closer. As a result, the LOOCV pixel has the weakest possible spatial correlation.

[0151] Experimentally, LOOCV MAE acts as an upper bound on LCP MAE, at least for low upsampling factors (2x and 4x) where local SiSIFUS performs decently. Local SiSIFUS MAE and LOOCV MAE may be positively correlated for the various window shapes (again, for low upsampling factors where local SiSIFUS works well), hence LOOCV MAE can be used to estimate the best window size to maximise the a-priori unknown local SiSIFUS MAE.

[0152] Figure 9 depicts an example computing system 1100. The following description of the computing system 1100 is provided by way of example only and is not intended to be limiting.

[0153] The example computing system 1100 includes a processor 1104 for executing software routines.

[0154] Although a single processor is shown for the sake of clarity, the computing system 1100 may also include a multi-processor system. The processor 1104 is connected to a communication infrastructure 1106 for communication with other components of the computing system 1100. The communication infrastructure 1106 may include, for example, a communications bus, cross-bar, or network.

[0155] The computing system 1100 further includes a main memory 1108, such as a random access memory (RAM), and a secondary memory 1110. The secondary memory 1110 may include, for example, a hard disk drive 1112 and / or a removable storage drive 1114, which may include a floppy disk drive, a magnetic tape drive, an optical disk drive, solid state storage or the like. The removable storage drive 1114 reads from and / or writes to a removable storage unit 1118 in a well-known manner. The removable storage unit 1118 may include a floppy disk, magnetic tape, optical disk, removable solid state storage (e.g. SD card) or the like, which is read by and written to by removable storage drive 1114. As will be appreciated by persons skilled in the relevant art(s), the removable storage unit 1118 includes a computer readable storage medium having stored therein computer executable program code instructions and / or data.

[0156] In an alternative implementation, the secondary memory 1110 may additionally or alternatively include other similar means for allowing computer programs or other instructions to be loaded into the computing system 1100. Such means can include, for example, a removable storage unit 1122 and an interface 1120. Examples of a removable storage unit 1122 and interface 1120 include a program cartridge and cartridge interface (such as that found in video game console devices), a removable memory chip (such as an EPROM or PROM) and associated socket, and other removable storage units 1122 and interfaces 1120 which allow software and data to be transferred from the removable storage unit 1122 to the computer system 1100.

[0157] The computing system 1100 also includes at least one communication interface 1124. The communication interface 1124 allows software and data to be transferred between computing system 1100 and external devices via a communication path 1126. In various embodiments, the communication interface 1124 permits data to be transferred between the computing system 1100 and a data communication network, such as a public data or private data communication network. The communication interface 1124 may be used to exchange data between a plurality of different computing systems 1100 that together form an interconnected computer network. Examples of a communication interface 1124 can include a modem, a network interface (such as an Ethernet card), a communication port, an antenna with associated circuitry and the like. The communication interface 1124 may be wired or may be wireless. Software and data transferred via the communication interface 1124 are in the form of signals which can be electronic, electromagnetic, optical or other signals capable of being received by communication interface 1124. These signals are provided to the communication interface via the communication path 1126.

[0158] As shown in Figure 9, the computing system 1100 further includes a display interface 1102 which performs operations for rendering images to an associated display 1130 and an audio interface 1132 for performing operations for playing audio content via associated speaker(s) 1134.

[0159] As used herein, the term "computer program product" may refer, in part, to removable storage unit 1118, removable storage unit 1122, a hard disk installed in hard disk drive 1112, or a carrier wave carrying software over communication path 1126 (wireless link or cable) to communication interface 1124. A computer readable medium can include magnetic media, optical media, or other recordable media, or media that transmits a carrier wave or other signal. These computer program products are devices for providing software to the computing system 1100.

[0160] The computer programs (also called computer program code) are stored in main memory 1108 and / or secondary memory 1110. Computer programs can also be received via the communication interface 1124. Such computer programs, when executed, enable the computing system 1100 to perform one or more features of embodiments discussed herein. In various embodiments, the computer programs, when executed, enable the processor 1104 to perform features of the above-described embodiments. Accordingly, such computer programs represent controllers of the computer system 1100.

[0161] Software may be stored in a computer program product and loaded into the computing system 1100 using the removable storage drive 1114, the hard disk drive 1112, or the interface 1120. Alternatively, the computer program product may be downloaded to the computer system 1100 over the communications path 1126. The software, when executed by the processor 1104, causes the computing system 1100 to perform functions of embodiments described herein.

[0162] It is to be understood that the embodiment of Figure 9 is presented merely by way of example. Therefore, in some embodiments one or more features of the computing system 1100 may be omitted. Also, in some embodiments, one or more features of the computing system 1100 may be combined together. Additionally, in some embodiments, one or more features of the computing system 1100 may be split into one or more component parts. It will be appreciated that the elements illustrated in Figure 9 function to provide means for performing the various functions and operations of the computer systems as described in the above embodiments. Although a few preferred embodiments have been shown and described, it will be appreciated by those skilled in the art that various changes and modifications might be made without departing from the scope of the invention, as defined in the appended claim Figure 9 depicts an example computing system 1100. The following description of the computing system 1100 is provided by way of example only and is not intended to be limiting.

[0163] The example computing system 1100 includes a processor 1104 for executing software routines.

[0164] Although a single processor is shown for the sake of clarity, the computing system 1100 may also include a multi-processor system. The processor 1104 is connected to a communication infrastructure 1106 for communication with other components of the computing system 1100. The communication infrastructure 1006 may include, for example, a communications bus, cross-bar, or network.

[0165] The computing system 1100 further includes a main memory 1108, such as a random access memory (RAM), and a secondary memory 1110. The secondary memory 1110 may include, for example, a hard disk drive 1112 and / or a removable storage drive 1114, which may include a floppy disk drive, a magnetic tape drive, an optical disk drive, solid state storage or the like. The removable storage drive 1114 reads from and / or writes to a removable storage unit 1118 in a well-known manner. The removable storage unit 1118 may include a floppy disk, magnetic tape, optical disk, removable solid state storage (e.g. SD card) or the like, which is read by and written to by removable storage drive 1114. As will be appreciated by persons skilled in the relevant art(s), the removable storage unit 1118 includes a computer readable storage medium having stored therein computer executable program code instructions and / or data.

[0166] In an alternative implementation, the secondary memory 1110 may additionally or alternatively include other similar means for allowing computer programs or other instructions to be loaded into the computing system 1100. Such means can include, for example, a removable storage unit 1122 and an interface 1120. Examples of a removable storage unit 1122 and interface 1120 include a program cartridge and cartridge interface (such as that found in video game console devices), a removable memory chip (such as an EPROM or PROM) and associated socket, and other removable storage units 1122 and interfaces 1120 which allow software and data to be transferred from the removable storage unit 1122 to the computer system 1100.

[0167] The computing system 1100 also includes at least one communication interface 1124. The communication interface 1124 allows software and data to be transferred between computing system 1100 and external devices via a communication path 1126. In various embodiments, the communication interface 1124 permits data to be transferred between the computing system 1100 and a data communication network, such as a public data or private data communication network. The communication interface 1124 may be used to exchange data between a plurality of different computing systems 1100 that together form an interconnected computer network. Examples of a communication interface 1124 can include a modem, a network interface (such as an Ethernet card), a communication port, an antenna with associated circuitry and the like. The communication interface 1124 may be wired or may be wireless. Software and data transferred via the communication interface 1124 are in the form of signals which can be electronic, electromagnetic, optical or other signals capable of being received by communication interface 1124. These signals are provided to the communication interface via the communication path 1126.

[0168] As shown in Figure 9, the computing system 1100 further includes a display interface 1102 which performs operations for rendering images to an associated display 1130 and an audio interface 1132 for performing operations for playing audio content via associated speaker(s) 1134.

[0169] As used herein, the term "computer program product" may refer, in part, to removable storage unit 1118, removable storage unit 1122, a hard disk installed in hard disk drive 1112, or a carrier wave carrying software over communication path 1126 (wireless link or cable) to communication interface 1124. A computer readable medium can include magnetic media, optical media, or other recordable media, or media that transmits a carrier wave or other signal. These computer program products are devices for providing software to the computing system 1100.

[0170] The computer programs (also called computer program code) are stored in main memory 1108 and / or secondary memory 1110. Computer programs can also be received via the communication interface 1124. Such computer programs, when executed, enable the computing system 1000 to perform one or more features of embodiments discussed herein. In various embodiments, the computer programs, when executed, enable the processor 1104 to perform features of the above-described embodiments.

[0171] Accordingly, such computer programs represent controllers of the computer system 1100.

[0172] Software may be stored in a computer program product and loaded into the computing system 1100 using the removable storage drive 1114, the hard disk drive 1112, or the interface 1120. Alternatively, the computer program product may be downloaded to the computer system 1100 over the communications path 1126. The software, when executed by the processor 1104, causes the computing system 1100 to perform functions of embodiments described herein.

[0173] It is to be understood that the embodiment of Figure 9 is presented merely by way of example. Therefore, in some embodiments one or more features of the computing system 1100 may be omitted. Also, in some embodiments, one or more features of the computing system 1100 may be combined together.

[0174] Additionally, in some embodiments, one or more features of the computing system 1100 may be split into one or more component parts.

[0175] It will be appreciated that the elements illustrated in Figure 9 function to provide means for performing the various functions and operations of the computer systems as described in the above embodiments. Although a few preferred embodiments have been shown and described, it will be appreciated by those skilled in the art that various changes and modifications might be made without departing from the scope of the invention, as defined in the appended claims.

[0176] The features disclosed in the foregoing description, or in the following claims, or in the accompanying drawings, expressed in their specific forms or in terms of a means for performing the disclosed function, or a method or process for obtaining the disclosed results, as appropriate, may, separately, or in any combination of such features, be utilised for realising the invention in diverse forms thereof. While the invention has been described in conjunction with the exemplary embodiments described above, many equivalent modifications and variations will be apparent to those skilled in the art when given this disclosure. Accordingly, the exemplary embodiments of the invention set forth above are considered to be illustrative and not limiting. Various changes to the described embodiments may be made without departing from the spirit and scope of the invention.

[0177] For the avoidance of any doubt, any theoretical explanations provided herein are provided for the purposes of improving the understanding of a reader. The inventors do not wish to be bound by any of these theoretical explanations.

[0178] Any section headings used herein are for organizational purposes only and are not to be construed as limiting the subject matter described.

[0179] Throughout this specification, including the claims which follow, unless the context requires otherwise, the word “comprise” and “include”, and variations such as “comprises”, “comprising”, and “including” will be understood to imply the inclusion of a stated integer or step or group of integers or steps but not the exclusion of any other integer or step or group of integers or steps.

[0180] It must be noted that, as used in the specification and the appended claims, the singular forms “a,” “an,” and “the” include plural referents unless the context clearly dictates otherwise. Ranges may be expressed herein as from “about” one particular value, and / or to “about” another particular value. When such a range is expressed, another embodiment includes from the one particular value and / or to the other particular value. Similarly, when values are expressed as approximations, by the use of the antecedent “about,” it will be understood that the particular value forms another embodiment. The term “about” in relation to a numerical value is optional and means for example + / - 10%.

[0181] References:

[0182] [1] Matplotlib, Interpolations for imshow.

[0183] [2] J. R. Lakowicz, Time-Domain Lifetime Measurements, in Principles of Fluorescence Spectroscopy, edited by J. R. Lakowicz (Springer US, Boston, MA, 2006) pp. 97-155.

[0184] [3] A. Chambolle, An algorithm for total variation minimization and applications, Journal of Mathematical imaging and vision 20, 89 (2004).

[0185] [4] R. K. Henderson, N. Johnston, F. Mattioli Della Rocca, H. Chen, D. Day-Uei Li, G. Hungerford, R. Hirsch, D. Mcloskey, P. Yip, and D. J. S. Birch, A 192x128 time correlated spad image sensor in 40-nm emos technology, IEEE Journal of Solid-State Circuits 54, 1907 (2019).

[0186] [5] R. Zhang, P. Isola, A. A. Efros, E. Shechtman, and O. Wang, The unreasonable effectiveness of deep features as a perceptual metric, in Proceedings of the IEEE conference on computer vision and pattern recognition (2018) pp. 586-595. [6] D. P. Kingma and J. Ba, Adam: A method for stochastic optimization, arXiv preprint arXiv:1412.6980 (2014).

[0187] [7] S. Fadnavis, Image interpolation techniques in digital image processing: an overview, International Journal of Engineering Research and Applications 4, 70 (2014).

Claims

Claims:1 . A computer-implemented method for upscaling the resolution of an image, comprising: obtaining a co-registered pair of images of different modalities, wherein the pair of images comprises a first image having a first resolution and a second image having a second resolution that is lower than the first resolution; determining a plurality of local correlation functions, wherein each local correlation function is associated with a respective sub-region of the second image and represents a correlation between pixel data of the respective sub-region with corresponding pixel data of the first image; and generating an upscaled version of the second image using the plurality of local correlation functions and the pixel data of the first image.

2. The computer implemented method of claim 1 , wherein each respective sub-region has a size of 5x5 pixels of the second image.

3. A computer-implemented method for upscaling the resolution of an image comprising: obtaining a co-registered pair of images of different modalities, wherein the pair of images comprises a first image having a first resolution and a second image having a second resolution that is lower than the first resolution; producing training data by: obtaining a first plurality of patches in the first image; and, for each patch: obtaining a pixel in the second image corresponding to the centre of the patch; training an algorithm using the training data, the trained algorithm being configured to reproduce or predict pixel-data of the pixel in the second image corresponding to the centre of any patch of the first image; and generating an upscaled version of the second image based on the trained algorithm and pixeldata from the first image.

4. A computer-implemented method for upscaling the resolution of an image, comprising: obtaining a co-registered pair of images of different modalities, wherein the pair of images comprises a first image having a first resolution and a second image having a second resolution that is lower than the first resolution; determining a plurality of local correlation functions, wherein each local correlation function is associated with a respective sub-region of the second image and represents a correlation between pixel data of the respective sub-region with corresponding pixel data of the first image; producing training data by: obtaining a first plurality of patches in the first image; and, for each patch: obtaining a pixel in the second image corresponding to the centre of the patch; training an algorithm using the training data, the trained algorithm being configured to reproduce or predict pixel-data of the pixel in the second image corresponding to the centre of any patch of the first image; andgenerating an upscaled version of the second image using: the plurality of local correlation functions; the trained algorithm; and the pixel data of the first image.

5. The computer implemented method of claim 3 or claim 4, wherein the algorithm is a deep neural network.

6. The computer-implemented method of any one of claims 3-5, wherein producing the training data further comprises rotating and / or reflecting the plurality of patches.

7. The computer-implemented method of any one of claims 3-6, wherein a size of the patches is pre-determined based on features in the first image.

8. The computer implemented method of any preceding claim, further comprising determining a predicted error of the upscaling method.

9. The computer implemented method of claim 8, further comprising comparing the error with a threshold, and generating an alert based on the predicted error.

10. The computer implemented method of any preceding claim, wherein generating an upscaled version of the second image comprises; constructing an intermediate image; optimizing the intermediate image based on the intermediate image, the first image, and at least one of the local correlation functions or the trained algorithm by varying the intermediate image; and outputting an optimized intermediate image as the upscaled version of the second image.11 . The computer implemented method of claim 10, wherein on the intermediate image is optimised based also on the second image.

12. The computer implemented method of claim 10 or claim 11 , wherein on the intermediate image is optimised based also on the total variation of the intermediate image.

13. The computer-implemented method of any one of claims 10-12, wherein the intermediate image is optimised by minimising a cost function based on the first image and second image, and wherein at least one of the first image or the second image is weighted.

14. The computer-implemented method of any preceding claim, wherein obtaining a coregistered pair of images comprises: obtaining a pair of un-registered images; and co-registering the pair of un-registered images.

15. The computer-implemented method of any preceding claim, wherein the first image is an intensity image.

16. The computer-implemented method of any preceding claim, wherein the data comprises data sampled at different times; and wherein the second image is a lifetime image.

17. The computer-implemented method of any preceding claim, wherein obtaining a coregistered pair of images comprises: obtaining data that is not in the form of a single image; and computing, from the obtained data, the second image.

18. The computer-implemented method of claim 16 or claim 17, wherein the data is obtained using a single-photon avalanche diode array.

19. An imaging apparatus for capturing co-registered images and upscaling the resolution of an image comprising: an objective lens configured to direct light from a subject; a beam-splitter configured to receive light from the objective lens and direct it to: reference imaging equipment configured to generate a first image of the subject; low-resolution imaging equipment configured to generate a second image of the subject; and a processor configured to execute a computer-implemented method according to any preceding claim.

20. A computer system comprising a processor configured to perform a method according to any preceding claim .21 . A non-transitory computer-readable storage medium comprising instructions which, when executed by a computer, cause the computer to carry out a method according to any preceding claim.

Images