1 Introduction: basic requirements for a combination of improved resolution with pattern recognition

Microscope image analysis is used in many fields of technology, science and medicine. In all these applications, a common problem is represented by the low resolution of objects of interest, so that detection and classification are difficult and, because of correspondent poor resolution, analysis of such objects can be lead to misunderstanding. As a consequence, pattern recognition (PR) and identification of objects of interest is a first necessary step in order to improve the quality of image analysis, and the second step is to improve the resolution from low level to high level, or equivalently, multiple scale resolution.

PR applied to the imaging techniques is a term denoting supervision methods developed by the combination of machine learning, artificial intelligence and data mining. PR methods can be categorized according to the type of learning procedure used to generate output value from input images. For example, supervised learning procedure assumes that a set of training data has been previously provided [1]. Then, a further learning processing generates a model that attempts to meet two general requests: (1) perform as well as possible on the training data and (2) generalize as well as possible to new data or other input images. The other one category is represented by the unsupervised learning methods. This class of pattern recognition methods assumes training data that have not been previously labelled, and tries to find inherent patterns in the data that can be used to identify a correct output value for new data instances [2]. A combination of the two categories is also possible, mixing labelled and unlabelled data in suitable ways, depending on the data in input, and where the data to be labelled are the training data [3]. Most common pattern recognition algorithms are probabilistic in nature; in that, they use statistical inference to find label for a given instance. In many applications, many probabilistic algorithms output a list of n-best labels with associated probabilities, instead of simply best label. For example, PR can be used for classification, which attempts to assign each input value to one of the given set of classes; in this case, the n-best labels should be fairly small (a binary choice for n is the best, yes or not, 1 or 0, etc). All the microscopic techniques require highly sophisticated pattern recognition supervised methods, for example, PR is crucial for the segmentation procedure of components to be identified and classified in a large image. Since PR strategies are strongly dependent on the application, we focused the attention on methods connected to the identification or organs and components inside biological systems like animal cells. For example, systemic analysis of subcellular protein localization provides a key for understanding gene functions and physiological condition of the cells. However, recognition of cell images of subcellular structures highly depends on experience and becomes the rate-limiting step when classifying subcellular protein localization. Several research groups have extracted specific numerical features for the recognition of subcellular protein localization, but these recognition systems are restricted to images of single particular cell line acquired by one specific imaging system and not applied to recognize a range of cell image sources.

Our imaging technique is based on the scanning probe microscope [4], a technique that measures the cells touching them with a probe approaching a resolution on the nanometre scale. On the contrary, with respect to other microscope technique, like for example, the confocal microscopy, the interaction between the microscopic tip probe and the sample is limited to the cell cytoskeleton, i.e., the cell surface, so pattern recognition methods must be addressed to identify skeleton components, like microtubules, microfilaments and intermediate filaments.

Once a skeleton components is identified, a procedure of super-resolution (enhanced resolution) could improve the data analysis and mining. The problem of enhancing resolution of a single image refers to the task of obtaining a high-resolution enlargement of a single low-resolution image, commonly recognized as single-frame super-resolution problem. This problem is intrinsically ill-posed because generally multiple frame images are combined to obtain final high-resolution images. Accordingly, to enhance resolution for a single image, strong prior information must be realized. This information is available in either the explicit form of a distribution defined on an image class or an implicit form of an example image which leads to example-based super-resolution image [5–7]. One diffuse approach for the super-resolution of low-resolution image has been characterized by algorithms including nearest neighbour (NN)-based estimation. In such approach, there are two distinct steps: in the first one, pairs of low-resolution and corresponding high-resolution image patches are collected. Then, in the second one, the super-resolution step, each patch of the given low resolution is contrasted to the stored low-resolution patches, and the high-resolution patch corresponding to the nearest low-resolution patch and satisfying a certain spatial neighbourhood similarity can be evidenced so that can be considered as the final output of the super-resolution image.

In this paper, we propose a simple and general-purpose method able to combine high-resolution procedure with the classification and identification of components of microscopy imaging. The method is the combination of two main procedures. In the first one, PR methods are applied to the image for promising components to be detected and small regions containing component of interest are extracted using a feature finder. In the second one stage, high-resolution algorithms are applied to such identified components in order to approach a multiple scales of resolution. However, the method is indeed to be applied to any microscopy technique, as well as any image recording device. We have focused the attention on biological systems, like animal cells, recorded with an atomic force microscopy; the application of our method to amphorae’s images recorded with a Microsoft Kinect is presented in supplementary material.

The paper is organized as follows. The next section is devoted to present the strategies and processing for enhancing the resolution in single-frame acquisition technique such as AFM. The subsequent section shows the pattern identification algorithm as used in this paper. In turn, the results obtained using the combination of PR and super-resolution methods are presented and discussed.

2 Enhanced resolution strategy in single-frame acquisition

The current work extends the framework of Kim and Kwon to obtain enhanced resolution images from a single-frame dataset of biological system images recorded with an AFM instrument [5, 6]. Here, we will detail with the procedure; the application on the AFM images will be discussed in the next section. The procedure for the super-resolution approach is as follows: first, the input from a single-frame image is interpolated into the desired scale using cubic spline interpolation. The interpolation generates an X-image where high-frequency details are still missed, and hence, they must be estimated based on its band frequency components (LX) extracted by applying a Laplacian to X-image. Then, the estimate (Y) can be added to X to produce the super-resolved image. A local patch-based regression is adopted for the estimation of LX. As the output patches overlap with their neighbours, a dataset of candidates for each new pixel location in a 3D image is generated, i.e., a Z image. Any candidate of such dataset is obtained based on many local realizations of input images in a final possible output, and it contains various partial information of the possible final high-resolution image. A single image with enhanced resolution is finally obtained as a convex combination for each pixel of the set of candidate pixels based on their estimated likelihood. To improve the visual quality, the results are post-processed based on an appropriate estimated prior [6].

The regression-based method is better than the interpolation despite many information acquired from the prior knowledge, in particular on the discontinuity (highly changing gradients) of the input image. Then, the prior information can be increased and stored for the regression methods using a segmentation of the input image, once patter recognition methods are employed to identify components or subcomponents of objects in the input, in our case, filaments and microtubule on the cell cytoskeleton.

During the regression processing, patch pairs are randomly sampled from a set of low-resolution and corresponding high-resolution approached patches. For a given set of training data points \(\{(\mathbf{x}_{1}, \mathbf{y}_{1}),{\ldots },(\mathbf{x}_{l},\mathbf{y}_{l})\}\subset \mathfrak {R}^{\mathrm{N}}\times \mathfrak {R}^{\mathrm{N}}\), we minimize the following regularized cost functional E for the regressor \(\mathbf{f}=\{f^{1},{\ldots },f^{N}\}\):

$$\begin{aligned} E\left( \mathbf{f} \right) =\frac{1}{2}\sum _{i=1,\ldots ,N} {\left( {\sum _{j=1,\ldots ,l} {\left( {f^{i}\left( {\mathbf{x}_j } \right) -y_j^i } \right) ^{2}+\lambda ^{i}\left\| {f^{i}} \right\| ^{2}} } \right) }\nonumber \\ \end{aligned}$$
(1)

where \(y_j =\left[ {y_{j,\ldots ,}^1 y_j^N } \right] ^{T}\) and denotes a \( {L}^{2}\) Hilbert space. The regressor f can be expanded using a class of function as

$$\begin{aligned} f^{i}\left( \theta \right) =\sum _{j=1,\ldots ,l} {a_j^i } k^{i}\left( {\mathbf{x}_{j} ,\theta } \right) \quad \hbox {for} \quad i=1,{\ldots },N \end{aligned}$$
(2)

and

$$\begin{aligned} \left\| {f^{i}} \right\| ^{2}=\sum _{m,n=1,\ldots ,l} {a_{m}^{i} } a_{n}^{i} k^{i}\left( {\mathbf{x}_{m} ,\mathbf{x}_{n} } \right) \quad \hbox {for} \quad i=1,{\ldots },N\nonumber \\ \end{aligned}$$
(3)

where \(k^{i}\) is a Gaussian kernel [5]

$$\begin{aligned} k^{i}\left( {\mathbf{x},\mathbf{y}} \right) =\exp \left( {-\frac{\left\| {\mathbf{x}-\mathbf{y}} \right\| ^{2}}{\sigma _k^i }} \right) \end{aligned}$$
(4)

Equation (1) represents a sum of individual convex cost functional for each scalar-valued regressor \(f^{i}\) and, in some cases, the minimum can be obtained analytically. An analytical solution requires the construction and the inversion of N kernel matrices like \(( {[ {\mathbf{K}_{( {m,n})}^i } ]_{l,l} =k^{i}( {\mathbf{x}_m ,\mathbf{x}_n }), \ \hbox {for} \ i=1,\ldots ,N} )\) during the training phase and \(N\times 1\) kernel evaluations in testing, making the problem prohibitive by a computational point of view. As demonstrated by Kim and Kwon, this problem can be managed looking at the optimality of the solution by typing the regularization parameter and the kernel parameter for the regressors (i.e. fixing the best choices for \(\lambda ^{i}\) and \(\sigma _{k}^{i})\) and finding the minimizer of cost functional (1) only within the span of a sparse basis set \(\{k(\mathbf{b}_{l},\cdot ),{\ldots }, k(\mathbf{b}_{s},\cdot )\}\) (with \(s<<l\)):

$$\begin{aligned} f^{i}\left( \cdot \right) =\sum _{j=1,\ldots s} {a_j^i } k\left( {\mathbf{b}_j ,\cdot } \right) \quad \hbox {for} \quad i=1,{\ldots },N \end{aligned}$$
(5)

Inserting the new set of functions (5) inside Eq. (1), the cost functional (1) is a differentiable function of the basis points \(\{\mathbf{b}_{1},{\ldots },\mathbf{b}_{\mathrm{s}}\}\) and a gradient-based optimization can be applied [5]. However, due to the nonconvexity of cost functional (1) with respect to the basis points \(\{\mathbf{b}_{1},{\ldots },\mathbf{b}_{\mathrm{s}}\}\), a direct gradient descent can meet local minima, so a good heuristic is required to initialize the solution. Kernel ridge regression can be utilized and the sparse basis set can be constructed based on the combination of the kernel matching pursuit and gradient descent, thereby reducing the time complexity of training and testing for regression. The problem of ringing artefacts can be resolved by exploiting a prior over image structure which takes into account the discontinuity of pixels value across edges to this problem is dedicated the next section. Here, we wish to emphasize that as the result of the patch-based regression step, approximately N candidates are generated for each pixel location by sharing the hyperparameters and basis points. The information gained by such N candidates in term of magnification factor is approximately 2, constructing an image by choosing the best and the worst candidates (best and worst in term of the distance to the ground truth that is coinciding with the prior, i.e., the best PSNR-supposed value) from each spatial location of Z.

The final estimation of the pixel value for an image location (x, y) is obtained as the convex combination of candidates given in the form:

$$\begin{aligned} I(x,y)=\sum _{i=1,\ldots ,N} {w_i } \left( {x,y} \right) Z\left( {x,y,i} \right) \end{aligned}$$
(6)

where

$$\begin{aligned} w_i (x,y)={\exp \left( {-\frac{d_i \left( {x,y} \right) }{\sigma ^{i}}} \right) }\bigg /{\sum _{j=1,\ldots ,N} {\exp } }\left( {-\frac{d_j \left( {x,y} \right) }{\sigma ^{j}}} \right) \end{aligned}$$
(7)

and

$$\begin{aligned} d_i \left( {x,y}\right) =\left| {PZ\left( {x,y}\right) ^{T}W_i}\right| \quad \hbox {for}\quad i=1,{\ldots },N \end{aligned}$$
(8)

where PZ(x, y) is a vector constructed by concatenating all columns of a spatial patch (of size \(R\times R\times N)\) of Z centred at (x, y), while the set of parameters \(\{W_{i}\}\) is optimized as a function of the patch-based regression results on Z for a subset of training images. Details on such optimization can be found in the section describing the results.

In such method, combining kernel ridge regression with kernel matching pursuit and gradient descent, the parameters to be tuned and optimized are few, essentially, the parameters \(\lambda \) and \(\sigma \), the number of basis s, the input patch R, for images \(256\times 256\) pixels or \(512\times 512\) as used in our paper; the values were generally \(R=64 (8 \times 8), N=36 \,\, (6\times 6), s=200\) or 300. The remaining hyperparameters are chosen based on error rates of super-resolution results for a set of validation images. The values are determined by trading the quality of the super-resolution images, and it has been generally observed that the performance of enhanced resolution increases as s increases.

3 Pattern identification algorithm

The classical approach of PR for the identification of objects inside images is based on textures representation, feature extraction and classification stage. A colour image, or a grayscale one, contains a finite set of uniform and texture domains. Such domains are composed from pixels with closed values, because represent a special case of uniformity. For example, in a colour image, a connected set of pixels with the same RGB value can represent a uniform region represented in a common colour. As a consequence first our approach with PR application to biosystem images, the preprocessing approach is made by image smoothing, edge enhancement or contrast adjustment. After such operations, we should be able to solve the recognition of a textured region by looking for the uniformity of such regions.

Once a textured region is identified, we have to find all its occurrences, which means the clustering of all image regions in a proper number of categories, depending on their uniformity or texture, so different textured regions can be grouped in a number of classes depending on the specific PR requirements. Generally, the identification of textured region is strictly connected to the segmentation problem. The image segmentation task consists of dividing the input image in a number of different objects, such that all the pixels from a segment have a common property based on the similarity criterion [7]. There are various segmentation methods, such as boundary-based approaches, which use edge detection for segments extraction, or amplitude thresholding approaches. Edges correspond to sharp variations of image intensity and convey vitally important information in an image. A variety of algorithms exist for edge characterization and detection such as the statistical methods, [8], difference methods [9], and curve fitting [10].

We have used the edge-detection method [11], and the edge detection is driven by the gradient methods, i.e., considering for the uniform regions a single characteristic value in the features vector. A region \(R (M \times N)\) of a given image I (with \(R\subset I)\) is uniform if the difference between the greatest pixel value and the smallest pixel value is below a small enough chosen threshold, i.e., \(\forall p \in R, {\vert } p\hbox {-mean}(R){\vert }<T\), where p is a pixel and T the threshold value. The value of T is chosen with a gradient-based value, and for any pixel of the preprocessed pattern (similar region) to be identified, we have computed the mean of its \(n \times n\) neighbourhoods and consider the obtained value as a feature vector. The block dimension n must have an odd value such that the current pixel can be its centre. In addition, a greater and a smaller value for n depending on the image noise’s amount and the uniform regions’ sizes. The value of n is set interactively by the user, but it should be chosen such that all the regions of interest be preserved (a large n could delete the small uniform regions). The patterns with their features values are then passed to the classification stage. Feature detection and extraction must be using the feature vector; for a R(x, y) image region with \(M \times N\), we define a moment m as follows

$$\begin{aligned} m_{i,j} =\sum _{i=1}^N {\sum _{j=1}^M {x^{i}y^{j}} } R\left( {x,y} \right) \,\, \hbox {with} \,\, i,j \,\, \hbox {integers} \end{aligned}$$
(9)

Once featured vectors have been identified, pattern classification can be done in a supervised manner or in a nonsupervised one. A supervised learning can be applied to a classification problem if the number of classes is known. If no information about classes is available or only a little is known about them, then the problem can be managed only with nonsupervised learning.

In a supervised classification approach, we must introduce a vector classifier (VC) and apply this classifier as a variant of the k-means algorithm for it [7, 12]. Let \(C_{1}, C_{2},{\ldots }C_{K}\) be the classes, K being the only information we know about them, and \(x_{1},{\ldots }x_{N}\), be the objects we want to cluster (the \(x_{i}\) correspond to the image pixels). Their feature vectors, \(V(x_{i}), i=1,{\ldots }, N\), are sorted in ascending order. Let \(\alpha _{1}=\hbox {min}V(x_{i}), \alpha _{K}=\hbox {max}V(x_{i}), i=1,{\ldots }, N\) and \(d( \alpha _{1}, \alpha _{K})\) the distance between them. Dividing this distance to K-1, we obtain K equally spaced vectors \(\alpha _{1}, \alpha _{2}, {\ldots }, \alpha _{K}, \hbox {d}(\alpha _{1}, \alpha _{K})/K-1\) being the distance between 2 successive vectors. For each i, min \(d(\alpha _{i}, V(x_{j}))\) is computed and the minimum Euclidean distance obtained is included in the training set. The feature vectors \(V(x_{i})\) are then substituted by a training feature vectors composed by K elements. The advantage of such method is that can be applied jointly to the image segmentation case, but also in a general pattern classification problem. The number of classes K is set interactively by displaying the image on the screen and the k-means algorithm computes successfully the \(C_{1}, C_{2},{\ldots },C_{K}\) classes, being classes of pixels, until an optimization criterion is reached. Optimization criterion is based on the rule that all the feature vectors are in the right place (moving from the centroids of the \(C_{i}\)) and no more position changes are needed.

Our pattern recognition strategies were applied to animal cell images recorded with an AFM. On the contrary, to optical microscope, an AFM takes an image of the sample interacting with the surface by means of a probe tip, and the interaction between the AFM probe tip and cell is conditioned by the elastic response of the sample. As a consequence, the AFM shows the cytoskeleton structure of a cell making possible to apply only supervised pattern recognition methods because all the cytoskeleton components such as microtubules, intermediate filaments, microfilaments are well known and, generally, rather easily identifiable.

In the next section, we present and comment the results obtained by applying the combination of PR and super-resolution methods to the identification of cytoskeleton components as displayed by AFM imaging.

Fig. 1
figure 1

On the left, source image: AFM image of an osteogenic hSMC, number of pixel 256\(\times \)256, real dimensions \(xyz \,\,\, 5\,\upmu \hbox {m}\times 5\,\upmu \hbox {m}\times 0.8 \upmu \hbox {m}\). The network of microtubules and filaments with dimensions falling in the range 100–1000 nm are evident. The red objects in the middle image summarize the application of criterion (13) for a PR identification of microtubules. On the right, the same microtubules after threshold segmentation (colour figure online)

Full size image

4 Combination of enhanced resolution and pattern identification: the case of human mesenchymal stem cells (hMSCs) imaged with AFM

In this section, we condensate the combined methods on pattern recognition and enhanced resolution to identify skeletal components and organs of hMSCs cytoskeleton as measured using an AFM. AFM is usually the most used SPM technique to investigate cell skeletons [13–15]. The AFM works using a probe to image the cell sample [4]. Such tiny probe can be considered as a paraboloid with a final sphere (the radius of the sphere falls in the range 10–20 nm) in permanent or intermittent contact (correspondent to two different modes of image acquisition) with the sample that can be consideredgenerally flat. Based on the contact force between the probe and the cell sample, the image recorded with the AFM presents a shot of the cell cytoskeleton, Fig. 1. Being the cytoskeleton composed by a complex network of different cell components, principally, microfilaments, intermediate filaments, microtubules, it can be a rather complex challenge to identify the different cyto-components.

In our procedure, the first step is represented by the application of PR approach to identify and extract components (possibly classified as microtubule, intermediate filaments or microfilaments). The goal is to identify a specific area \(I_{0}(x, y)\subseteq I(x,y)\) of the source image I, corresponding in our case to a microtubule. The PR algorithm based on the edge-detection method operates as follows:

  1. 1.

    in the selected area \(I_{0}\) identifies a possible a central pixel \(p_{0}\) that can be considered as a barycentre pixel of \(I_{0}\);

  2. 2.

    then, it applies a kind of region growing method based on the gradient value of pixel neighbourhood.

  3. 3.

    All the pixels in a neighbourhood region are correlated, and a local estimate of the correlation surface is made.

In particular, the third step of out edge-detection approach has the advantage for a fast computation since the value of the gradient itself address us where to look forward. Assuming \(I_{0}(x,y)\) as the target image of size \(I_{a} \times I_{b}\), and the W (x, y) is a windowing function (of the same size of \(I_{0}\) or less) containing the object of interest, whose gradient goes to zero at the edges (a sort of boundary condition), then the normalized squared correlation is given by (for sake of simplicity in continuous form)

$$\begin{aligned} C^{2}\left( {\rho ,\theta } \right) =\frac{\left[ {{\int } {\int {I\left( {x,y} \right) \left( {I_0 W} \right) \left( {x-\rho ,y-\theta } \right) \mathrm{d}x\mathrm{d}y} } } \right] ^{2}}{{\int } {\int {I^{2}\left( {x,y} \right) W\left( {x-\rho ,y-\theta } \right) } } } \end{aligned}$$
(10)

with the assumption that \({\int } {\int {\left( {I^{2}W} \right) \left( {x,y} \right) } } \mathrm{d}x\mathrm{d}y=1\). Taking the derivative of the normalized correlation we find

$$\begin{aligned} \nabla C^{2}\left( {\rho ,\theta } \right)= & {} -2\frac{\left[ {{\int } {\int {I\left( {x,y} \right) \left( {I_0 W} \right) \left( {x-\rho ,y-\theta } \right) \mathrm{d}x\mathrm{d}y} } } \right] \times \left[ {\int {\int {I\left( {x,y} \right) \nabla \left( {I_0 W} \right) \left( {x-\rho ,y-\theta } \right) \mathrm{d}x\mathrm{d}y} } } \right] }{{\int } {\int {I^{2}\left( {x,y} \right) W\left( {x-\rho ,y-\theta } \right) } } } \nonumber \\&\quad +\left[ {\frac{{\int } {\int {I\left( {x,y} \right) \left( {I_0 W} \right) \left( {x-\rho ,y-\theta } \right) \mathrm{d}x\mathrm{d}y} } }{{\int } {\int {I^{2}\left( {x,y} \right) W\left( {x-\rho ,y-\theta } \right) } } }} \right] ^{2}\times \left[ {{\int } {\int {I\left( {x,y} \right) \nabla W\left( {x-\rho ,y-\theta } \right) \mathrm{d}x\mathrm{d}y} } } \right] \end{aligned}$$
(11)

that, at the least, leaves four terms to be calculated:

$$\begin{aligned}&{\int } {\int {I\left( {x,y} \right) \left( {I_0 W} \right) \left( {x-\rho ,y-\theta } \right) \mathrm{d}x\mathrm{d}y} } \nonumber \\&{\int } {\int {I\left( {x,y} \right) \nabla \left( {I_0 W} \right) \left( {x-\rho ,y-\theta } \right) \mathrm{d}x\mathrm{d}y} } \nonumber \\&{\int } {\int {I^{2}\left( {x,y} \right) W\left( {x-\rho ,y-\theta } \right) } } \nonumber \\&{\int } {\int {I\left( {x,y} \right) \nabla W\left( {x-\rho ,y-\theta } \right) \mathrm{d}x\mathrm{d}y} } \end{aligned}$$
(12)

These equations can be easily transformed to the discrete domain, reducing the calculation of the correlation gradient at a central pixel p and the correspondent 8 neighbours. In addition, a further simplification can be made considering that \(\nabla W\) and \(\nabla I_{0} W\) can be precomputed when the shape and dimension of the target are previously well known (highest information available by the prior), or a reasonable threshold value is taken in consideration, before to starting the search [16]. Once the PR procedure is applied, the algorithm generates a connected region corresponding to a microtubule by adding a pixel connected with a surrounding pixel selected with the criterion of the derivative of the correlation function with a threshold value. In practice, starting from the selected pixel \(p_{0}\in I_{0}\), a new pixel in \(I_{0}\) if and only if:

$$\begin{aligned}&p_{i+1} \in I_0 \Leftrightarrow p_i \in I_0 ,\quad \left| {p_{i+1} -p_i } \right|<\Delta _{R}\nonumber \\&\left| {p_{i+1} -p_0 } \right| <\Delta _{A} \end{aligned}$$
(13)

One main problem applying the criterion (13) is the identification of overlapping objects. To discriminate two overlapping objects, we apply a supplemental identification in z-direction, this can be done analogously to (13) using two values in \(\hbox {z}_{\pm }\), one for the top-object, (i.e. \(p_{i+1} \in I_0 \Leftrightarrow p_i \in I_0 ,\left| {p_{i+1} -p_i } \right| >z_+ )\), and the other one for the bottom object (i.e. \(p_{i+1} \in I_0 \Leftrightarrow p_i \in I_0 ,\left| {p_{i+1} -p_i } \right| <z_- )\). The results of such supplemental criterion along with the criterion (13) are shown in Fig. 1.

The enhancing resolution procedure is adapted by Kim and Kwon, as described before, and it consists of the following steps: first, we interpolate the inputs into the appropriate scale using cubic spline interpolation. Second, the high-frequency details which are missing during the interpolation processing are estimated based on the band frequency components extracted by applying the Laplacian to the interpolated data.

One basic problem on the treatment of the image is the presence of sharp changes that can be identified as the edges. The edges can be smoothed, or oscillation occurs to compensate applying some averaged value to descent gradient methods as given by the criterion (13) when we calculate Eq. (11). Generally, the success of such procedure needs suitable average values and such choice requires the combination of enhanced or super-resolution methods with PR classification. This requirement is satisfied by an ad hoc algorithm recently introduced by us: the PRIAR algorithm [17].

The PRIAR (pattern recognition image augmented resolution) is a tool we have developed that gives the possibility to enhance the resolution of low-resolution images. It enhances the spatial resolution and the depth of the image. The implemented algorithm follows the scheme: (1) image classification; (2) blind super-resolution on single frame; (3) pattern-analysis; (4) reconstruction of the discovered pattern.

PRIAR runs on MATLAB 2013b, and the main features of the computer used to test the program as follows: CPU Intel(R) Core(TM) i5-2400 CPU @ 3.10 GHz (4 cores) with a L2 cache of 6144 KB with 16 GB of RAM. The program has been developed as mono-task, the idea is to have a pipe-flow where each piece calculates a subproblem sequentially. PRIAR works with as input a low-resolution image and produces as output a high-resolved computed image where the objects to be identified are recognized and substituted with their correspondent models. The execution time depends on the size of the input data. A typical execution time for patch images of \(128\times 128\) pixels, 8 bit per pixel, is of about 70 s. The occupied memory depends on the image size, and the presented typical example needs of about 256 MB of RAM.

The \(\hbox {MATLAB}\circledR \) code that describes algorithm implementation is detailed elsewhere [17]; here, we resume the basic functions included inside with the following pseudo-code:

figure a

The algorithm models a general-purpose method to analyse the images; in fact, it uses a modular solution that permits to extend the class of the recognized images and the image enhancement.

The line 02 of the pseudo-code calls a data-mining function that classifies the input image. In this particular implementation, the function distinguishes only two classes of images: images representing a grating and images representing cells. The classification is important during the next stages to choose the better segmentation to adopt.

The line 07 reports the function call that super-resolve the image. During this step, the analysed image (that is a low-resolution image) is improved by increasing the colorimetric and spatial resolution. The PRIAR implements the following SR algorithms: Kim and Kwon [5], spline interpolation [18], nearest-neighbour interpolation [19], bilinear interpolation [20], bicubic interpolation [21], box-shaped kernel interpolation [22], Lanczos-2 kernel interpolation [21] and Lanczos-3 kernel interpolation [22]. We made a comparison with different enhancement methods. The test images come from the BSDS300 dataset of Berkeley [23]. The experimental results show that our method enhances the texture regions, while it maintains the contrast of image. The results are freely available on our webpage [24]. A comparison between our SR method and Kim-Kwon method is presented in supplementary material.

The line 08 describes the first step to map the object of interest (microtubule). In each case, it needs a seed. The seed can be simply a point of the object we are going to map or a polyline that follows a part of the object we are going to discover.

The main idea of the explore algorithm is to crawl the surface until certain parameters are respected. The role of the parameters can be summarized, in fact, that the gap between two next pixel must be enough small. At the end of this exploration, we have an initial path that is a subset of the object we are going to trace. The hypothesis that the initial path is a subset of the final path is guarantee by the strong constraint that are applied by the explore algorithm. The complexity of this algorithm is linear with the number of the pixel that composed the matrix. The explore algorithm performs a local exploration so it is a local search algorithm.

Line 09 shows the edge discover function call. The function used with a simple pattern, like a grating image, is different by the one adopted with cell images, complex pattern, in order to avoid unwanted shadows. In fact, in order to edge the grating super-resolved image, we combined the Sobel method with the Laplacian of Gauss method [1, 25]. In the case of complex pattern like cells, we combined the Canny method with the Zero-Cross method [1, 25]. To demonstrate the performance of the proposed method, we compare our approach with Chan-Vese segmentation method (error % of matching edge patterns) on reconstructed patterns using a conventional peak SNR (PSNR) measurement [26, 27], Table 1.

Table 1 Comparison between PRIAR and Chan–Vese method based on the determination of Error % with different PSNR
Full size table
Fig. 2
figure 2

Segmented microtubules as in Fig. 1, but with the application of PRIAR processing for enhancing resolution with a magnification of a factor 2. On the left: 3D representation of the tubules as input image. On the right: 3D representation of the correspondent enhanced resolution image with edges treatment as given by Eq. (11) and denoising processing

Full size image

Line 10 describes the final step that recognize the object to be tracked. It combines the initial path with the segments obtained by the edging procedure. In our model, we are searching something that has the shape of a linear uniform pipe displaced over a ground so the function calculates some characteristic of our initial path (e.g. the mean high, the orientation and its centroid) and tries to match them with the characteristics of the segmented areas. When they match, a new set of pixel is joined to the initial one. At the end of this function, the object is tracked and the function provides a path that maps the discovered object. It is important notice that in contrast to the explore algorithm, this is a global search algorithm.

Line 11 shows the application of the last function that enhances the image using the model of the previously identified object that coincides with the prior. In particular, we underline that the lines 1–10 use a standard approach to image analysis that consists of improving resolution of the input image and identify a pattern over it. The line 11 changes the algorithm prospective because the recognized pattern is used to improve the image reversing the roles of PR and super-resolution.

In order to define a pattern including all the possible object classifiable as microtubules, the algorithm is applied again on the surrounding pixels of the object already identified as a microtubule using Eq. (11).

Once the objects of interest are identified, their image quality is enhanced of a magnification factor 2, they can be managed to reproduce them in a 3D rendering visualization for image analysis. This rendering generally requires a supplemental processing involving edge-detection procedure and possible spikes eliminated with smoothing as in Fig. 2.

5 Conclusive remarks

In this paper, we have described a simple and general-purpose method able to combine high-resolution procedure with the classification and identification of patterns from cells images recorded with an atomic force microscopy. The method is composed of two different stages combining pattern recognition as a preparatory step to the enhanced resolution of the identified components. First (pattern recognition), promising components (objects of interest) are detected in the input image and small regions containing such objects of interest are extracted using a feature finder. Second, high-resolution algorithms are applied to such identified components in order to approach a multiple scales of resolution. In turn, the method has been built and indeed to be applied to any microscopy technique.