## Materials and Methods
### Methodology  
As mentioned above, the classification of multi-temporal images is affected by many factors, such as missing data in cloudy images, optimum training sample size, and appropriate temporal resolution for image acquisition that limits the accuracy. Moreover, extraction of appropriate features (e.g., the temporal gradient) from the temporal signatures of crop types, along with suitable times such that differentiate crops are two major issues that should be considered in the multi-temporal image classification. In this regard, different features, like VIs, the temporal gradient of VIs, and spectral bands, are taken into account as new features. To address the aforementioned issues, we developed the following methodology (Figure 4). **Step 1:** A set of images with an amount of cloudiness less than 20% in their metadata, were selected from input time-series images. Radiometric correction is then applied to the selected images (i.e., 10 images for this study area) and cloudy pixels in each image are set to no data as missing data for all bands. Hence, for each pixel, a feature vector with a dimension of 80 elements (8 \(\times\) 10) is achieved. **Step 2:** As can be seen in Figure 4, the missing data retrieval is a preprocessing task before RF classification. Indeed, the missing data of pixels due to cloudy periods are very important, since many pixels are ignored during multi-temporal image classification. The proposed method for missing data retrieval is described in detail in Figure 5. \begin{table}
\begin{tabular}{c c c c c} \hline \hline
**Class** & **\# of Training Fields** & **\# of Test Fields** & **\# of Training Pixels** & **\# of Test Pixels** \\ \hline
**Wheat** & 128 & 280 & 559 & 1265 \\
**Barley** & 15 & 51 & 93 & 307 \\
**Other** & 121 & 252 & 1419 & 3027 \\ \hline \hline \end{tabular}
\end{table} 
Table 2: The number of training and test data in the study area. Figure 3: General agronomic calendars of the crops in the study area. **Step 3:** After missing data retrieval, the VIs are generated. Previous studies have not investigated the impact of VIs on classification of wheat, but it is clear that RF was most effective at detecting them. In the current study nine indices, are computed based on the corresponding equations listed in Table 3, as new features. The goal of considering 9 VIs is to create features that enhance classification accuracy. To retrieve the missing data, all the bands of a pixel corresponding to cloudy periods in a given training/test data set (i.e., wheat) that contain no data were removed from the feature vector (i.e., an image which captured in time 3) and flagged as cloudy pixels. The corresponding bands from the feature vector of the cloudless pixels of that class were also removed (i.e., the corresponding bands of cloudy times are removed in cloudless pixels of that class) to obtain the same dimensionality with the cloudy pixels (i.e., 72 \(\times\) 1). In the reduced multi-dimensional feature space, for each cloudy pixel (i.e., the green pixel in Figure 5), k-nearest cloudless pixels (i.e., the red pixels in Figure 5) are then found. Therefore, the bands of pixels with no data were reconstructed using the mean of the k-nearest neighbor cloudless pixels, to obtain a feature vector with full dimensionality. In other words, the missing data of cloudy pixel bands are filled with the mean of the k-nearest cloudless pixels, which are spectrally similar to cloudy pixels in the feature space. This was also applied for all of the cloudy pixels of the scene after the appropriate multi-temporal classification model was selected. Figure 4: Flow chart describing the method developed in this study for multi-temporal image classification. **Step 4:** Since the temporal signature plays an important role in the interpretation of the growing steps of crops, a temporal gradient is used to quantitatively describe the seasonal change.