Thread Counting in Plain Weave for Old Paintings Using Regression Deep Learning Models

International Journal of Computer Vision https://doi.org/10.1007/s11263-025-02473-9 Thread Counting in Plain Weave for Old Paintings Using Regression Deep Learning Models Antonio Delgado1·Juan José Murillo-Fuentes1·Laura Alba-Carcelén2 Received: 22 September 2023 / Accepted: 3 May 2025 © The Author(s) 2025 Abstract In this paper, we introduce a novel algorithm designed to improve thread density estimation in canvas analysis. Our approach incorporates three major contributions. First, we eliminate the need for post-segmentation processing by integrating regression techniques, enabling the deep learning (DL) model to directly compute thread density. This does not only reduce computational time but also shifts the training focus from locating crossing points to minimizing thread counting errors, thereby enhancing accuracy. We develop and rigorously evaluate various models, selecting the one with optimal performance through a hyperparameter search. Second, we refine the data generation process by dynamically adjusting filter lengths based on initial thread density estimates and incorporating equalization. We also enhance data augmentation. Third, we implement semi-supervised training to expand the dataset and fine-tune model weights. This involves incorporating new inputs into the training set when both the DL model and Fourier transform yield similar density estimates for new paintings. Our proposed algorithm demonstrates superior performance in thread density error reduction and operational efficiency compared to previous DL segmentation solutions for masterpieces from Ribera, Velázquez, or Poussin. Additionally, it has been effectively applied to identify fabric matches between canvases attributed to different authors, showcasing its practical applicability in art analysis. Keywords X-ray image processing ·Canvas weave ·Thread Counting ·Deep Learning ·Regression ·Inception 1 Introduction 1.1 Interest in Fabric analysis In the forensic study of a painting, the analysis of the fabric plays an important role, see Alba-Carcelén et al. (2020)for further details. Within the features of the fabrics, the number of threads per centimeter (cm) in both vertical and horizontal Communicated by Dimosthenis Karatzas. BJuan José Murillo-Fuentes [email protected] Antonio Delgado [email protected] Laura Alba-Carcelén [email protected] 1Departamento de Teoría de la Señal y Comunicaciones, ETSi Universidad de Sevilla, Camino de los Descubrimientos sn, 41092 Sevilla, Spain 2Gabinete de Documentación Técnica, Museo Nacional del Prado, C. de Ruiz de Alarcón, 23, 28014 Madrid, Spain orientation (density of threads) is usually studied, especially in plain weave. Plain weave fabrics are the support for a vast number of paintings thanks to their good compromise between robustness and simplicity. This type of cloth is characterized by an orthogonal intertwining of the warp and weft threads, see Fig. 1.a. The separation between the warp threads follows a deterministic pattern that is given by their placement on the loom while the weft usually follows a random arrangement. In fabric analysis with image processing and computer vision, the X-ray of the canvas is usually used. In Fig. 1.b we include some examples of 1 cm side patches of X-ray plates from different canvases with plain weave. Note that in these patches we cannot identify warp and weft. Accordingly, we will simply refer to the vertical and horizontal threads, as seen in the image. Nowadays, with the introduction of automatic thread counting (Johnson et al., 2010,2013), the curators not only avoid the tedious task of counting but they have access to a map of densities for vertical threads throughout the painting and another for horizontal ones. Because the separation of the threads in plain weave is not regular and depends on the 123 International Journal of Computer Vision Fig. 1 In (a) a sketch of plain weave with warp and weft yarns in weaving as described in Barlow (1878) and in (b) patches of 1 cm side from X-ray plates of paintings using plain weave manufacturing and the loom, if we find the same shared pattern in the maps of two canvases, we can conclude that both come from the same bolt. Therefore, by matching thread density maps curators can accomplish precise studies about the timing or authorship of masterpieces. There are many examples of the application of thread counting to forensic studies. See for example the Automated Canvas Analysis (2019) where the paintings by Van Gogh are analyzed and sequenced. Van Gogh ordered a new canvas only when the previous roll was completely used. Since the type of fabric differs from one bolt to another, by matching the canvas we can sequence them and at the same time corroborate authorship. A similar project is being developed for the work by Vermeer, see Johnson (2017). So far, we have reported applications of thread counting to study if some paintings share the same fabric. Thread counting is also applied to check for integrity, i.e., if some original parts of the canvas are missing or to analyze transformations of the size of the paintings through time. 1.2 Motivation and Contributions Given the patches from different fabrics included in Fig. 1.b, it is evidenced that while counting threads may be an easy task in some cases, it can be really hard in others. Poor resolution, noise from the painting itself, cracks, distortions such as rotations, or a high density of prime make it difficult to count the number of threads in many scenarios. We have recently proposed a method for thread density estimation based on deep learning (DL) (Goodfellow et al., 2016) and in particular on the U-Net model (Ronneberger et al., 2015), by segmenting crossing points in Delgado-Bejarano et al. (2022a,b,2023). Segmentation of crossing points based on DL presented good results compared to the state-of-the-art, i.e., Fourier transform (FT) based tools. However, it focused on the location of the crossing points between threads. The estimation of the densities from the segmentation result is computed in posterior processing. Furthermore, in the segmentation DL paradigm, it is quite difficult to link the segmentation result with the error in thread counting to train the DL network. One of the major difficulties of the DL is the need for labeled samples in the training stage. This training is run once and then the weights of the model are fixed and ready to be used in the analysis of any painting. A rich and wide set of samples, covering different qualities of fabrics, ranges of thread densities, and noises are needed. However, the labeling stage is not only time-consuming but quite hard to accomplish for high-density thread fabrics and noisy images. To create a large dataset, data augmentation (DA) (Shorten & Khoshgoftaar, 2019) was exploited. You might see DelgadoBejarano et al. (2023) for further information on DA, where we cropped areas of large labeled samples to generate inputs. Besides, we further augmented the data set by randomly rotating the result of the cropping. Also, the images were pre-processed. A filter of fixed length was used in this step. We propose a novel algorithm that incorporates three major contributions. First, we avoid the post-processing of the segmentation result by resorting to regression and forcing the DL model to directly compute the thread density itself. The advantages are twofold. On the one hand, we remove the processing after the segmentation, reducing computational running time. On the other, we use as a training loss function the error in the estimation of the thread density, i.e., the model is not trained to locate crossing points, but to minimize the error in the thread counting, providing a more accurate result. We develop and analyze different models to retain the one with the best performance where we perform an optimized search as described in Snoek et al. (2012); O’Malley et al. (2019) to set their hyperparameters. As a second enhancement, we deeply review the data generation process. In the pre-processing step, a linear-regression-based algorithm has been developed to automatically adapt the pre-processing parameters for each X-ray plate based on a rough estimation of its thread density. Also, histogram equalization is incorporated into the pre-processing stage to enhance thread highlighting and facilitate counting. In the DA, we limit the random rotation applied to some of the images and increase the number of images by 40% by also cropping central parts of the labeled samples. Third, we introduce semi-supervised training to increase the dataset and further refine the weights of the model: when processing a new full painting, inputs for which both the DL and the FT provide similar density estimates are incorporated into the training dataset. As a result, the proposed approach outperforms the previous approaches in terms of thread density error and reduces the running time compared to the DL segmentation solution. The organization of the paper is as follows. In Section 2we described related previous works. The methodology is developed in Section 3. Improvements in the DA and pre123 International Journal of Computer Vision processing are developed in Section 4. The proposed models, based on DL regression, are included in Section 5. These models are trained and their performance is compared in Section 6. A semi-supervised model is presented in Section 7. In Section 8we evidence the good performance of the proposal compared to previous approaches for several analyses. Section 9is devoted to conclusions. 2 Related works In the literature, we find three different main techniques for thread density estimation. Namely, those based on frequency analysis, feature extraction followed by machine learning, and DL. 2.1 Frequency Analysis In Escofet et al. (2001); Johnson et al. (2013); Simois and Murillo-Fuentes (2018); Murillo-Fuentes and Alba (2018) the authors describe the analysis of fabrics based on the discrete 2D-FT. In Johnson et al. (2010,2013), this idea was exploited to obtain thread counting maps, by applying the discrete 2D-FT to patches of images all over the X-ray plate and finding the maxima. In Fig. 2we include a brief description of the process. For a given area in the original image, the 2D-FT is computed. Then the position of the maximum in the horizontal axis provides the vertical thread frequency estimation. The value is depicted in the corresponding position in the vertical density map, where each value indexes a different color. The main advantage of the discrete 2D-FT approach is that it is an unsupervised method. Besides, it is usually robust to the presence of noise, such as the painting itself, and artifacts in the X-ray such as nails or stretchers. In scenarios where we have a quite uniform warp-weft pattern of threads, the FT provides quite a good result. However, in Delgado-Bejarano et al. (2023) it is reported that the FT fails whenever we have different distances between nearby threads, the width of the threads varies, or the warp is tighter than the weft. In Fig. 2we observe that the final result is quite noisy, as the vertical thread density map should exhibit well-defined vertical lines. You may compare this result to the one in Fig. 5, obtained by using the model later proposed in this work. 2.1.1 Feature Extraction Approach Another approach that has been proposed to study the thread density in canvases is based on feature extraction and machine learning (Maaten & Erdmann, 2015). The approach was presented as an automatic thread-level canvas analysis tool, hereafter denoted by ATCA. It extracts histograms of oriented gradient features. These features are the input to Fig. 2 Description of the method to compute the vertical thread density map with the FT. Examples of the outputs of the blocks are included below arrows machine learning methods such as support vector machines and Bayesian logistic regression classifiers to determine if a location in the image is a crossing point or not. By applying the approach through every pixel, it is possible to generate thread density maps. This method exhibits a major drawback: it is necessary to label a large number of crossing points in the X-ray plate image of the to-be-analyzed canvas. This is cumbersome for the practitioner. Besides, the method depends on a set of parameters whose optimal values depend on the image itself. While some of them are automatically computed, others might need manual adjustment. 2.1.2 Deep Learning Algorithms based on DL have been applied to several problems in art, e.g. Sizyakin et al. (2020); Roberto et al. (2020); Pu et al. (2020); Zou et al. (2021) where convolutional neural networks (CNN), U-Net (Ronneberger et al., 2015), and auto-encoders (AE) (Rumelhart et al., 1986) were exploited. Within the scope of this work, in Delgado-Bejarano et al. (2022a,b,2023) the authors proposed to segment crossing points using DL, to later estimate the thread densities. In Fig. 3we describe the process. A DL-based segmentation model locates crossing points ina1cmsidesquare patch from the X-ray image. Then image processing is used to estimate the average distance between crossing points, i.e., the thread density. The process is repeated for locations along the canvas to get the density maps. The pre-processing involves double filtering of fixed length kernels to reduce the difference in the mean of the input images and increase the contrast. Note that the labeling is performed just for the training set. Hence, compared to the ATCA method, this method avoids any labeling of the to-be-analyzed painting. The algorithm provided good estimates where FT methods failed. In the models in Fig. 2to Fig. 5, the definition of the pre-processing module is open. In the frequency analysis, there is no relevant literature on this issue. While in Delgado-Bejarano et al. (2023) a pre-processing approach was proposed. In this manuscript, we further develop this one. 123 International Journal of Computer Vision Fig. 3 Description of the method to compute the vertical thread density map with DL-based segmentation. Examples of the outputs of the blocks are included below 3 Methodology The pipeline of the novel proposed solution has four main steps. The first one is the generation of the dataset. In Fig. 4we include a description of the data generation process. Every Xray plate available is pre-processed. Pre-processing involves a resampling to 200 pixels per centimeter (PPCM), if needed, and some filterings to enhance the image. Then we perform a random cropping to get samples of size 300×300. Each sample is annotated. With DA we obtain several instances from each sample, with the corresponding annotations, of sizes 200×200. Finally, from the annotation, we get the estimation of the densities, e.g., for the vertical threads. We can use the spatial counting approach developed in Delgado-Bejarano et al. (2023) to provide this density estimation. This approach detects adjacent crossing points to estimate horizontal and vertical distances between them. It can be applied to the annotations of the crossing points in the labeled image samples. The resulting instances are stored in the dataset. In the second stage of the pipeline, we define our model. A DL model with convolutional layers followed by a fully connected (FC) network is proposed, that directly provides the estimation of the vertical thread counting for an input image. The third step in the pipeline involves training. We present instances from the dataset to the model, compute the loss function using the estimation and the label, and update the weights of the network using gradient descent. In the last stage of the pipeline, we face how to analyze an X-ray plate image. The image is divided into patches that, after pre-processing, are fed sequentially to the model. The outputs with the vertical threads density estimations for every patch are arranged into a matrix, where each entry corresponds to the location of the input patch. The matrix is usually represented as an image, where estimates for threads countings are translated into colors for a better perception. In Fig. 5we include a workflow of the process. At this point, it is interesting to underline that compared to the solution in Fig. 3the deep learning model directly estimates the thread density from the input. We do not need to perform the density estimation step. For the same reason, in the learning stage, the loss function directly accounts for Fig. 4 Data Generation Process: 1) the X-ray plate image is preprocessed, 2) then cropped to get some samples, of 1.5 cm side, 3) each sample image is labeled to mark the position of the threads, 4) for every sample image and its annotation, we generate several pair of instances for the dataset, of 1 cm side Fig. 5 Flowchart of the DL-based Regression approach for thread counting. Examples of the outputs of the blocks are included below the arrows the error in the thread counting, instead of the error in the segmentation. Besides, note that once we estimate the vertical thread densities, by repeating the procedure with the inputs rotated by 90◦we have the estimation of the horizontal ones. Accordingly, in the training dataset, we include all instances and their 90◦rotated counterparts, with the corresponding horizontal thread countings as labels. In the fourth step, in production, we consider further refining the weights of the model in a semi-supervised approach. Then, after re-training to refine the weights of the model, we finally estimate the output. In the remainder of this work, we describe these steps in detail, including some examples of application to some case studies. 4 Data Generation We used the same X-ray plates and random samples as in Delgado-Bejarano et al. (2023): 240 labeled samples of 1.5 cm side from 36 selected paintings from the Museo Nacional del Prado (MNP). Canvases from Rubens, Velázquez, Lorena, Swanevelt, Dughet, Poussin, Both, Lemaire, and Ribera, among others, were included in the dataset. The fabrics of these paintings have several densities in the range of 6 to 23 threads per cm (thr/cm), different resolutions of the image, and several noise conditions. This 123 International Journal of Computer Vision encompasses the usual densities found in canvases, see for example the analysis of thread densities in Carbonnel (1980) for French painting. The data generation follows the procedure described in Fig. 4. We next describe the pre-processing and the DA blocks. 4.1 Pre-Processing We have adapted the pre-processing in Delgado-Bejarano et al. (2023). We first resample the input image to have 200 PPCM if needed. Then perform a mean and local standard deviation filtering, see Delgado-Bejarano et al. (2023)fora full description. Finally, the pixel values were adjusted to cover the whole available range. We have incorporated two major improvements to this procedure: 1) we propose an automatic and adaptive estimation of the dimensions of the filter, i.e., the kernel, for each X-ray plate and 2) we introduce a histogram equalization to emphasize the threads and to facilitate thread counting by the neural network. In the mean and standard deviation filterings we use as a filter a square kernel of size kto estimate the mean or the standard deviation around a pixel to update its value with the result. In Delgado-Bejarano et al. (2023) both filters had fixed sizes. We found that images with a better contrast could be obtained by adapting the kernel sizes to the mean fabric thread density. In practice, we want the kernel to have a width large enough to cover an area to include all ranges of values and at the same time small enough to adapt to changes within the image. For this reason, we propose to apply a dynamic window based on the average distance between threads of each fabric. First, we need to estimate thread density values throughout the analyzed painting. Second, we provide an expression to estimate the filter size, based on the largest distances found between threads in pixels. The method to estimate the size of the kernel, k,topreprocess the image is described in Alg. 1. An average thread density value can be quickly extracted with a coarse frequency analysis. We first pre-process the X-ray image with an initial value, k0=21. We used the 512 points FT of 1 × 1 cm patches from the X-ray plate sampled every 7-15 cm, depending on the canvas size. Once the thread densities are estimated for this set of samples, we define the parameter t by averaging the largest values above the mode for the vertical and horizontal density histograms. Then, we divide the value by 1.5, to improve the estimation in cases of extreme counts. We used linear regression to provide the optimal kernel size, that approximates the averaged distance between threads, in pixels. The coefficients of the linear regression model in (1) were adjusted by using a large set of pairs of count values (t) - kernel size (k) obtained from previous frequency studies. For low values of t, e.g. 6 threads per cm, we have k=31.65 in (1), which is approximately the distance Fig. 6 Input image obtained with (a) fixed window size and (b) dynamic window size. We have highlighted some areas and one of them has been zoomed, where after using adaptive kernel size the contrast of the image has been improved between threads, in pixels. For high values, e.g. 25, we have k=14.55. The method is run twice, to refine the estimation. If k1is the value obtained for kin the first iteration, we assume that k0=k1in the second iteration. We assume this value to be good enough to provide an accurate thread density estimation with FT for the second iteration. Algorithm 1 Kernel Size Estimation Input: Scanned X-ray plate image and initial kernel width k0=21. 1: Pre-process the X-ray with a kernel width of k0. 2: Run the FT approach to estimate the horizontal and vertical thread densities, for some set of croppings of the X-ray, in vectors thand tv, respectively. 3: Compute the mode: ˜ th=mode(th) ˜ tv=mode(tv). 4: Compute the average above the mode: th=mean(th(i):th(i)>˜ th) tv=mean(tv(i):tv(i)>˜ tv) and t=((th+tv)/2)/1.5. 5: Estimate k, k=−0.90 ·t+37.05.(1) 6: Repeat Steps 1-5 with k0=k. 7: Round kto the nearest value to 14.55, and then to the lowest odd number. 8: Output: k In Fig. 6.a we include a pre-processed sample with fixed size kernel and in Fig. 6.b we have the same sample preprocessed with adaptive kernel size. We highlighted three areas where a better contrast is observed, i.e., it is easier to distinguish crossing points from dark areas between threads. We also zoomed in on one of them. The other enhancement we have introduced is related to the transformation of the pixel values at the end of the 123 International Journal of Computer Vision pre-processing. We propose an off-the-shelf equalization approach as follows. We estimate the histogram, h(x),for the gray levels of the input image, in unsigned integer 8 bits format, normalize it to sum 255, and compute its integral, h(x), and use the result as a look-up table. For a given input value, x, the output yields y=h(x). This is applied as a final step before entering the DL regression model. 4.2 DA As described in Sec. 3, we used the same X-ray plates and random samples as in Delgado-Bejarano et al. (2023), which are 240 labeled samples of 300x300 pixels and 200 PPCM, i.e., sample images of 1.5 cm side. In Delgado-Bejarano et al. (2023), we defined a DA to exploit the labeled samples and obtain a much larger dataset. In this process, we cropped different 200x200 pixels (1 cm side) areas and transformed them to create several instances from every sample as follows. First, four instances were generated by cropping the four corners of the image. Then, four rotated versions of two patches given by the (50:250, 50:250) and (65:265, 65:265) pixels of the sample were added, with random angles in the ranges [2, 7], [-2, -7], [8, 12], [-8, -2]. Two additional rotated versions of patches of the sample images were added, given by the (65:265, 65:265) pixels, with random angles in the ranges [2, 7], [-2, -7]. The dataset was further augmented by repeating the whole process after a) left-right flipping, b) updown flipping, c) rotating 90◦, d) rotating 90◦plus left-right flipping, and e) rotating 90◦plus up-down flipping the full sample. As a result, there are 60 instances of size 200x200 for every labeled 300x300 pixel sample. Since we labeled 240, we obtained 14400 images. We generated twice the instances for some paintings, to enhance the balance of the dataset. Paintings with extreme, low or high, thread densities were selected. By repeating the generation of samples for these fabrics, we increased the number of images to 21540. In this work, two improvements have been proposed to the DA in Delgado-Bejarano et al. (2023). The first one is related to the rotations. The fabric of the canvases is usually distorted in some areas. In fact, along the stretchers, the fabric bends due to the tightness around the nails. In other areas of the canvas, the threads are not perfectly arranged vertically or horizontally; at some points, the deviations can be severe. In Delgado-Bejarano et al. (2023) randomly rotated images were included in the dataset to allow the model to learn the densities in rotated scenarios. However, we found that the maximum allowed rotations were larger than needed, and the model parameters might be biased towards angle deviations not found in practice. Accordingly, we have limited the maximum random rotations applied to half the maximum in Delgado-Bejarano et al. (2023), and random rotations applied in the DA now vary in the range [-6◦,6 ◦]. We also investigated enhancing the DA by increasing the number of instances per every 1.5 x 1.5 cm labeled sample. In Delgado-Bejarano et al. (2023), the cropped instances were taken from areas in the corners of the labeled samples. We have added the cropping of the central regions given by the (50:250, 50:250), (65:265, 65:265), (35:235, 35:235), and (80:280, 80:280) pixels and the flipped versions and 90◦rotated version, i.e., +24 instances per labeled sample. No randomly rotated versions of these +24 images are included. In summary, we now get 84 instances of 1x1cmfor every labeled sample, 40% more instances than in DelgadoBejarano et al. (2023), 240 x 84 = 20160 instances overall. By including the extra instances to balance the dataset, the final number of instances rises to 29820. These instances are divided into training (23949), validation (4368), and test sets (1512). All instances coming from one painting are assigned to the same subset. The test subset is not used in the training stage or to select the model or hyperparameters but to analyze the final results. 5 Proposed Model We designed some models to perform regression, providing the vertical density of threads per cm. In this section we present the one with best performace, used later to analyze some canvases. In the following, models are described by using a flowchart where arrows indicate a layer and rectangular prisms represent the intermediate features. In the convolutional layers, unless otherwise indicated, we use a basic layer denoted hereafter as the inception block, see Fig. 7. This block follows the inception paradigm in DelgadoBejarano et al. (2023) as it was shown to be useful when dealing with images of different densities of threads. By applying the inception paradigm, we have convolutional kernels of three sizes in parallel within the same layer: 3 ×3, 5×5, and 7 ×7, see Fig. 7. The results of the convolutions are concatenated at the output. Therefore, the number of features at the output is the number of filters multiplied by 3. In the diagrams of the models, a horizontal arrow denotes a layer of convolutional inception blocks. The number above the arrow is the number of filters used. The number below the prisms are of the form width ×height ×features.After convolution, within the inception block, we apply batch normalization followed by a ReLU activation function. In the diagrams, we depict a horizontal dashed line to indicate that the upper part is a CNN-based model while the lower one is a dense network. Vertical arrows in the CNN part indicate a 2D max-pooling, while in the dense network, they are FC layers. In between the CNN and dense networks, we have a flattening to adapt the output of the CNN network to the input of the dense part. FC layers use ReLU as activation functions except for the output neuron that uses a linear activation. 123 International Journal of Computer Vision Fig. 7 Inception block: The input is convolved with nisets of filters of size 3 ×3, 5 ×5, and 7 ×7. The results of these ni×3 convolutions are concatenated at the output Fig. 8 Inception VGG-based regression model 5.1 Inception VGG-Based Regression Model We propose a model for regression borrowing from the VGG16 (Simonyan & Zisserman, 2015) architecture. Three inception blocks are used in the deeper convolutional layer. The number of dense layers and the number of neurons in each of them are defined following the VGG architecture. There are two dense layers with 64 times the number of filters of the first layer, 8 ·64, plus one output layer of size one. An optimized search (Snoek et al., 2012; O’Malley et al., 2019) has been used to set the optimal number of convolutional filters to 8. The full model has over 9 million parameters. 6 Training and Testing Results 6.1 Other Proposed Models We proposed other designs, as follows, that will be tested in this section. We first designed a model from scratch, based on an architecture with an encoder followed by an Fig. 9 Inception regression model (Reg) FC. Therefore, this architecture has two main parts. The first part performs the feature extraction via convolutional layers. Each convolutional layer has two consecutive inception blocks and 2D max-pooling is applied between layers to reduce dimensionality. The second part uses a flattening layer and a dense network that performs the regression according to the extracted features. We performed an optimized search (Snoek et al., 2012; O’Malley et al., 2019) to set the optimal value for the number of convolutional and dense layers, the number of filters in each layer, and the number of neurons in each dense layer. The optimal architecture has an encoder with 6 double CNN layers, i.e., each layer has 2 consecutive inception blocks. At the output of the encoder, the model includes 5 dense layers with 100, 100, 80, 100, and 1 neuron. The number of kernels at each CNN layer is described in Fig. 9, it doubles as we go down in the encoder, except for the last convolutional layer. The full model has over 10 million parameters. This model will be referred to as Reg. We also proposed a design corresponding to the one in Fig. 9but in every inception block, we add the input to the output. The full model has over 10 million parameters. This model will be denoted by Res. For the sake of comparison, we also trained the Segmentation DL model in Delgado-Bejarano et al. (2023), including the improvements in Section 4. This method will be denoted by Seg. Finally, we also developed a model based on the Segmentation DL. This model follows the inception regression architecture (Reg), but the convolutional part was taken from the Segmentation DL (Delgado-Bejarano et al., 2023), with its weight values initialized to the ones of the Segmentation DL. It will be referred to as Reg-Unet. 123 International Journal of Computer Vision 6.2 Setup The models described above were programmed using Python 3.9, CUDA 11.2, and Keras-Tensorflow 2.5.0. The input grayscale image size was 200 ×200 ×1 and the default batch size was 32. We used the Adam optimizer (Kingma & Ba, 2015) with a default learning rate value of 1·10−3.Early stopping is used with a latency of 65 epochs and a limit of 450 epochs. The code was run on an NVIDIA RTX A6000. We trained every model 10 times. In every training, the weights of the models were randomly initialized. The datasets were also different throughout the training, as in the DA different random rotations were applied. The same datasets were used in all models. The normalized mean absolute error (NMAE) was used as the loss function, where the error was normalized by the value of the label, i.e., the thread density: L(y,ˆ y)=1 m m  i=1 |yi−ˆyi|/yi(2) where mis the number of instances per batch, yis a vector with the labels for the instances, and ˆ yis a vector with the corresponding estimations. In each training, we randomly initialized the weights of the model and randomly set the rotations in the DA to generate the training and validation datasets (see Delgado-Bejarano et al. (2023) for further details). 6.3 Training and Testing In Fig. 10 we include the box diagram with the results of the validation loss for the regression approaches developed: the VGG-based regression model in Section 5.1 denoted by R-VGG, the inception regression model Reg) , and its residual version (Res), the adaptation of the segmentation model (Reg-Unet), and the segmentation based approach in Delgado-Bejarano et al. (2023) (IncDice). Given the validation NMAE values in Fig. 10, it can be concluded that the Reg-Unet approach is not a good option, i.e., reusing the encoder of the segmentation is not useful. On the contrary, the Reg, R-VGG, and Reg-Res approaches provide good low values with low variation of the results for the 10 trainings. The difference in the mean and variances of the validation loss are not relevant for these three models, Friedman’s test provides a value of p=0.202, indicating that we cannot conclude significant differences between the validation losses obtained. The Segmentation DL model in Delgado-Bejarano et al. (2023) with the lowest NMAE achieved 1.12% while with the improvements in Section 4,wehaveanNMAEof 1.04%, i.e., we have a 0.08% improvement, but significantly above the regression-based approaches. Fig. 10 Box diagram of validation loss, i.e., the average NMAE (%) of the density estimations for the validation set, the 10 trainings and the 4 regression models considered, and the Segmentation DL. Outliers are represented with diamonds Among the evaluated methods the R-VGG exhibits the lowest value of NMAE for the validation set, 0.92%. We will use this model in the analysis of the canvases later in this work. Hence, compared to the segmentation DL approach the regression DL method has a 0.12% lower error and does not need further signal processing as it directly computes the thread densities. When evaluated in the test set, see Fig. 11, regression DL models have a much better error value than the Segmentation DL. Specifically, the NMAE in the test has decreased from 1.51% of the Segmentation DL to 1.02% of the R-VGG for the best model in validation. Friedman’s test between R-VGG and Segmentation DL provides a value of p = 0.011. The NMAE of the Segmentation DL as in DelgadoBejarano et al. (2023), i.e., without the improvements in the generation of the dataset, was 1.61%. 6.4 Ablation Study To analyze the impact of improvements in Section 4we perform an ablation study where for the model with the lowest validation loss, the R-VGG, we train and report the value of the NMAE for the test set when each of the proposed improvements is not used. We include in Fig. 12 the NMAE for the R-VGG regression and this same model with no equalization (N-Eq), without taking central crops (C.Crops) in DA, when maximum random rotation angles in the DA are doubled (Angle), and using kernels of fixed size (FW). Finally, we include the NMAE when pre-processing is avoided and original instances are directly fed into the R-VGG model (Orig.) 123 International Journal of Computer Vision Fig. 11 Box diagram of the average NMAE (%) of the density estimations for the test set, the 10 trainings and the 4 regression models considered and the Segmentation DL. Outliers are represented with diamonds Fig. 12 Box diagram of the average NMAE (%) of the density estimations for the test dataset in the ablation study: the regression model considered (R-VGG) without equalization (N-Eq), a pre-processing with fixed window (FW), removing images in DA taken from the central parts of the labeled samples (C.Crops), limiting the maximum value of the random angle in DA (Angle) and no using pre-processing (Orig). In each case, 10 different trainings were run In Delgado-Bejarano et al. (2023) large and low values in input images were clipped and then re-scaled to exploit the full available range. Here we introduce histogram equalization that, in view of Fig. 12, helps to improve the estimation of the thread densities. In practice, we found that the improvement is noticeable for some X-ray plates but irrelevant for others. On the other hand, we also improve the NMAE using a variable window for the mean and standard deviation filtering in pre-processing. Note that we will have larger reductions in the error for fabrics with optimal window sizes different from the fixed value used in Delgado-Bejarano et al. (2023). Therefore, the variable window size is especially interesting for canvases with extreme values for thread density values. Regarding DA, and given Fig. 12, either using more samples (C.Crops) or limiting the maximum random rotation (Angle) does not greatly improve the results, however, a slight gain is obtained when focusing on the best NMAE value. Note that these techniques involve a residual percentage of the overall set of instances, hence the impact is bounded. Finally, it can be concluded given the NMAE of the model with no pre-processing but exploiting DA, see Orig in Fig. 12, that the performance very much deteriorates. Hence, preprocessing is quite needed. 6.5 Errors in Test Instances For the sake of completeness, in Fig. 13.a we included the result of the normalized absolute error (NAE) in the vertical () and horizontal (•) density estimation for 72 instances of the test set. These instances come from 3 different canvases. Instances number 0-23, 24-47, and 48-71 correspond to the first, second, and third paintings, respectively. We used the model with the lowest validation error, i.e., the R-VGG. In Tab. 1we report the lowest NMAE obtained in validation and test for the models In Fig. 13.a we included the errors for the regression DL R-VGG model (R-VGG). In Fig. 13.b and Fig. 13.c, we depicted the error for the Segmentation DL as in Delgado-Bejarano et al. (2023) and the FT, respectively. The error of the FT is quite large for the instances coming from the first and third canvases. Compared to the Segmentation DL model, the error of the regression DL decreases for almost all the test samples. It can be observed in Fig. 13.c that the FT fails to estimate horizontal thread densities for the first canvas, instances 0-23. It also poorly estimates vertical thread densities for instances 54, 55, 60-65, 68, and 69. For the other ones, it provides a good estimation although the average NMAE is above the DL-based approaches. The error is quite low for some instances, e.g. 28, 29, 35, 43, 49, 51, 56, 58, 66, 67, 70, and 71. In some of these, the NMAE is lower than in the DL-based approaches, e.g. 29, 35, 51, or 70. By comparing the Regression DL, in Fig. 13.a, and the Segmentation DL, in Fig. 13.b, we conclude that, in general, the NMAE of every instance is reduced. With the proposed Regression DL model we reduce from 22 to 9 the number of instances with NMAE above 2.5%. Besides, it significantly reduces the NMAE for the horizontal thread densities in instances 0-23. While the reduction in instances 24-47 is more evident in the vertical thread densities. 123 International Journal of Computer Vision Johnson, D.H., Jr., C.R.J., & Erdmann, R.G. (2013). Weave analysis of paintings on canvas from radiographs. Signal Processing 93, 527–540 Tilborgh-(contact), L. (2019). Automated Canvas Analysis: a new tool for sequencing Van Gogh’s paintings. https://www. vangoghmuseum.nl/en/about/knowledge-and-research/researchprojects/automated-canvas-analysis Johnson, C. Richard, J., Ed, W.A.S. (2017). COUNTING VERMEER: USING WEAVE MAPS TO STUDY VERMEER’S CANVASES. https://countingvermeer.rkdstudies.nl/contents/ Goodfellow, I., Bengio, Y., & Courville, A. (2016). Deep Learning. Cambridge, MA, USA: MIT Press. Ronneberger, O., Fischer, P., & Brox, T. (2015). U-Net: Convolutional networks for biomedical image segmentation. In: Lecture Notes in Computer Science, 9351.https://doi.org/10.1007/978-3-31924574-4_28 Delgado-Bejarano, A., Murillo-Fuentes, J.J., & Alba-Carcelén, L. (2022). Crossings Segmentation in Plain Weaves for X-rays of Canvases with Deep Learning. In: 7th IP4AI Meeting Computational Approaches for Technical Imaging in Cultural Heritage. National Gallery Delgado-Bejarano, A., Murillo-Fuentes, J.J., & Alba-Carcelén, L. (2022). Crossings Segmentation in Plain Weaves for X-rays of Canvases with Deep Learning: Technical Details. In: 7th IP4AI Meeting Computational Approaches for Technical Imaging in Cultural Heritage. National Gallery Delgado-Bejarano, A., Alba-Carcelén, L., & Murillo-Fuentes, J.J. (2023). Crossing Points Detection in Plain Weave for Old Paintings with Deep Learning. Accepted in Eng. Appl. Artif. Intell. https:// arxiv.org/abs/2302.11924 Shorten, C., & Khoshgoftaar, T. M. (2019). A survey on image data augmentation for deep learning. Journal of Big Data, 6(1), 1–48. Snoek, J., Larochelle, H., & Adams, R.P. (2012). Practical bayesian optimization of machine learning algorithms. In: Advances in Neural Information Processing Systems, 25 O’Malley, T., Bursztein, E., Long, J., Chollet, F., Jin, H., Invernizzi, L., et al. (2019). KerasTurner. https://github.com/keras-team/kerastuner Escofet, J., Millán, M. S., & Ralló, M. (2001). Modeling of woven fabric structures based on Fourier image analysis. Applied Optics, 40, 6170–6176. https://doi.org/10.1364/ao.40.006170 Simois, F. J., & Murillo-Fuentes, J. J. (2018). On the power spectral density applied to the analysis of old canvases. Signal Processing, 143, 253–268. https://doi.org/10.1016/j.sigpro.2017.08.006 Murillo-Fuentes, J.J., & Alba, L. (2018). In: Bastidas, D.M., Cano, E. (eds.) Thread Counting in X-Rays of Plain-Weave Painting Canvas, pp. 91–105. Springer, Cham. https://doi.org/10.1007/978-3-31975316-4_7 . Maaten, L., & Erdmann, R. G. (2015). Automatic thread-level canvas analysis: A machine-learning approach to analyzing the canvas of paintings. IEEE Signal Process. Mag., 32(4), 38–45. https://doi. org/10.1109/MSP.2015.2407091 Sizyakin, R., Cornelis, B., Meeus, L., Dubois, H., Martens, M., Voronin, V., & Pizurica, A. (2020). Crack detection in paintings using convolutional neural networks. IEEE Access, 8, 74535–74552. https:// doi.org/10.1109/ACCESS.2020.2988856 Roberto, J., Ortego, D., & Davis, B. (2020). Toward the automatic retrieval and annotation of outsider art images: A preliminary statement. In: 1st Int. Work. on Artificial Intelligence for Historical Image Enrichment and Access (AI4HI), Marseille, France Pu, W., Sober, B., Daly, N., Higgitt, C., & Daubechies, I., Rodrigues, M.R.D. (2020). A connected auto-encoders based approach for image separation with side information: With applications to art investigation. In: IEEE Int. Conf. on Acoustics, Speech and Signal Process. (ICASSP), pp. 2213–2217 Zou, Z., Zhao, P., & Zhao, X. (2021). Virtual restoration of the colored paintings on weathered beams in the forbidden city using multiple deep learning algorithms. Advanced Engineering Informatics, 50,. https://doi.org/10.1016/j.aei.2021.101421 Rumelhart, D., Hinton, G., & Williams, R. (1986). Learning internal representations by error propagation. Cambridge: MIT Press. Carbonnel, K. V. (1980). A study of french painting canvases. Journal of the American Institute for Conservation, 20(1), 3–20. Simonyan, K., & Zisserman, A. (2015). Very deep convolutional networks for large-scale image recognition. In: Bengio, Y., LeCun, Y. (eds.) 3rd International Conference on Learning Representations, ICLR 2015, San Diego, CA, USA, May 7-9, 2015, Conference Track Proceedings. http://arxiv.org/abs/1409.1556 Kingma, D.P., & Ba, J. (2015). Adam: A method for stochastic optimization. In: 3rd Int. Conf. on Learning Representations, ICLR 2015, May 7-9, San Diego, California, USA. http://arxiv.org/abs/ 1412.6980 Aradillas, J. C., Murillo-Fuentes, J. J., & Olmos, P. M. (2021). Boosting offline handwritten text recognition in historical documents with few labeled lines. IEEE Access, 9, 76674–76688. https://doi.org/ 10.1109/ACCESS.2021.3082689 Ribera, J. (1632). Ixion. Museo Nacional del Prado (P001114) Poussin, N. (1636). The Triumph of Pan. National Gallery (NG6477) Poussin, N. (1636). The Triumph of Silenus. National Gallery (NG6477) Maaten, L. (2015). Canvas Analysis. https://lvdmaaten.github.io/ canvas/ Velázquez, D.R. (1632). Antonia de Ipeñarrieta y Galdós and her Son, Luis. Museo Nacional del Prado (P001196) Velázquez, D.R.: (1632) Diego del Corral y Arellano. Museo Nacional del Prado (P001195) Rizi, F. (1660). An Artillery General. Museo Nacional del Prado (P001127) Mozo, F.D.H.E.: Pope St. Leo I the Great. Museo Nacional del Prado (P007113) (17th Century) Fondón-García, I., Simois, F.J., & Murillo-Fuentes, J.J. (2014). Software tool for thread counting in X-rays of plain-weave painting canvas. In: Int. Conf. on Non-Destructive Investigations and Microanalysis for the Diagnostics and Conservation of Cultural and Environmental Heritage (ART2014), Madrid, Spain Publisher’s Note Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations. 123