Phenology-Based Unsupervised Rapeseed Mapping Using Multitemporal Data

With a large field of view, remote sensing is a useful technique that can provide a new way to acquire crop acreage and spatial distribution. As an important oil crop, rapeseed is widely planted in many countries. Besides, rapeseed mapping is of great significance for food security and policy regulation. The objective of this work is to explore an automatic and effective framework for rapeseed mapping using multitemporal multispectral images (MSIs) captured by multiple sensors. The automatic rapeseed mapping framework consists of two stages. Concretely, color transfer is introduced and combined with vegetation indexes to detect rapeseed for the first time. After that, the most general supervised classification methods are applied to optimize the initial mapping result. In the experiments, the performance of different methods is evaluated on multisensor data within four years (2017–2020). Moreover, to analyze the effect of temporal information, the mapping performance of using different MSIs as input is systematically compared. Experimental results illustrate the effectiveness of the proposed framework, in which color transfer makes an important role. Some valuable findings are also obtained which can be further used for global rapeseed mapping.

grow [6]. In another word, the phenological information of crops can be recorded by remote sensing techniques.
Rapeseed is one of the major oil crops and widely planted in China, which can be used to produce edible oil and biofuels [7], [8]. Monitoring the temporal and spatial dynamics of rapeseed is meaningful for the rapeseed sector, national food security, and other relevant stakeholders [9]. However, the annual rapeseed production and yield data are acquired by artificial statistics (i.e., field surveys and producer reports) in many regions currently. This method is time consuming and laborious, which also cannot obtain the rapeseed distribution map with detailed spatial information. It is necessary to realize automatic rapeseed mapping using remote sensing images.
According to the previous studies, rapeseed mapping methods can be roughly divided into three strategies: special phenological-period-based methods; multitemporal-based methods; and spectrum-based methods. These methods are introduced as follows.
1) Special phenological-period-based methods: The special phenological-period-based methods leverage the crop growth stages that contain more discriminative information to identify rapeseed. At the flowering stage, the prominent yellow flower is an important feature of rapeseed, which can be used to distinguish rapeseed from other crops and yield prediction. For example, Li et al. [10] pointed out that the best period for rapeseed detection based on Landsat TM images is the flowering period in Shou County, Anhui Province, China. Wang et al. [11] mapped rapeseed at the provincial scale by extracting HSV (hue, saturation, and value) and spectral features from Chinese Gaofen satellite no. 1 (GF-1) multispectral images acquired at six continuous flowering stages. Moreover, Han et al. [12] first detected the flowering and pod phases of different regions, then utilized the spectrum at the flowering period and the scattering signal at the pod period to identify rapeseed with multiple thresholds. Meng et al. [13] attempted to select the optimal temporal window for wheat and rapeseed mapping in Zhongxiang city, Hubei Province, China. According to this research, Sentinel-2 (S2) images from the middle and later stages of the growth cycle are conducive to crop mapping .
Considering the flowering period is significant, some studies try to extract yellow flowers and monitor the flower data by designing a rapeseed index. For example, Sulik et al. [14] proposed the ratio yellowness index (RYI) to indicate the yellowness in a canopy. Later, Sulik et al. [7] designed the normalized difference yellowness index (NDYI), which is a better indicator of yield potential and monitoring the peak flowering period than This work is licensed under a Creative Commons Attribution 4.0 License. For more information, see https://creativecommons.org/licenses/by/4.0/ the normalized difference vegetation index (NDVI). To enhance the weak flower signal, Zang et al. [15] proposed the enhanced area yellowness index (EAYI) based on the peak of difference yellowness index (DYI) and the valley of NDVI during the estimated flowering period.
2) Multitemporal-based methods: Multitemporal-based methods utilize time-series remote sensing corresponding to the entire growth cycles of plants as input. Compared with the aforementioned method, this kind of method makes full use of phenological features and can extract more temporal spatial information for crop detection. Since time sequence data of Moderate Resolution Imaging Spectroradiometer (MODIS) is stable, lots of studies utilized MODIS products to identify crops. For example, Tao et al. [9] built a decision tree with MODIS time-series data to extract the rapeseed distributed in the Middle Reaches of the Yangtze River Valley in China. In addition, Meng et al. [5] also proved that using time-series images can obtain better rapeseed mapping results than only using images at the flowering stage. Besides, Ashourloo et al. [16] analyzed the spectral and temporal profile of various crops, then developed a canola index (CI) to separate rapeseed from other crops in the flowering stage.
3) Spectrum-based methods: It is acknowledged that hyperspectral image (HSI) is able to provide abundant spectral and spatial structure information [17], [18], [19]. The multirange spectral feature also can be used to identify rapeseed. For example, a spectral feature fitting method is performed on Hyperion data to extract rapeseed planting areas in [20]. Meng et al. [5] compared the performance of different data (HSI, monotemporal and multitemporal multispectral images) with different mapping methods, such as convolutional neural network (CNN), support vector machine (SVM), and random forest (RF). Based on this research, using HSIs as input can get a comparable result to using multitemporal and multispectral images.
Although many vegetation index (VI)-based and special phenological-period-based methods are developed, the performance of these methods depends on the imaging time and quality of flowering data to a certain extent. When the images are captured at the peak flowering stage, these methods work well. Unfortunately, due to the weather, revisit period of satellite and phenological differences, the image at the peak flowering stage is not always available. Moreover, they cannot be directly used for multiple years or over a large region because the stage of flowering differs temporally and spatially. In contrast, multitemporal-based methods focus on the temporal variation of crop canopy reflectance rather than the characteristics of the specific phenological stage. Since multitemporalbased methods make full use of phenological information, the crop mapping results are satisfactory when training samples are adequate. The spectrum-based methods can achieve a comparable result with abundant spectral information, but the imaging area of HSIs is limited. It is not suitable for large-scale crop identification. In summary, it is still a huge challenge to automatically identify rapeseed with a finer resolution over a large region where cloud-free high/medium-resolution data are limited.
To solve the sample size problem, this article proposed an unsupervised rapeseed mapping method that can detect rapeseed in the region where the croplands are distributed unevenly with small sizes and the cloud-free imagery data are restricted. In the study region, images with a relatively higher spatial resolution are appropriate for analyzing rapeseed. So, experimental data are derived from multisensor with medium spatial resolution, such as GF-1 wide-field view imager (WFV), Sentinel-2 Multi Spectral Instrument, Landsat 7 Enhanced Thematic Mapper Plus(L7 ETM+), Landsat 8 Operational Land Imager (L8 OLI), and Chinese Gaofen satellite no. 6 (GF-6) WFV. With the time-series MSIs, color transfer (CT) and VIs were combined to identify rapeseed roughly at first, which is an unsupervised process. Then, different classifiers that have been demonstrated effective in land-use and land-cover classification were introduced to optimize the initial rapeseed detection result. More in detail, the training and validation samples were randomly selected from the initial result to train the aforementioned classifiers. These classifiers were performed on the multitemporal MSIs subsequently. According to the experiment results, the unsupervised rapeseed mapping framework was analyzed from three aspects, including the effectiveness of the unsupervised rapeseed initial detection method with CT, the systematic comparison of the universal classifiers, and the influence of temporal information on the rapeseed mapping result. The contributions of this work are as follows.
1) A two-stage unsupervised rapeseed mapping framework is proposed, which consists of an initial detection and an optimization stage. The experiments are designed to validate the feasibility of this framework with multitemporal and multisensor MSIs from 2017 to 2020. 2) To the best of our knowledge, it is the first time to introduce color transfer for rapeseed mapping. Experiments illustrate color transfer can further enrich the initial detection result.
3) The influence of temporal information on rapeseed mapping is explored in this work. The amount of available cloud-free multitemporal MSIs is limited, it is necessary to find an appropriate combination of images as input. Thus, the relationship between the image of phenological periods and mapping results are studied. The remainder of this article is organized as follows. Section II introduces the proposed framework. Experimental settings and details of the study area are reviewed in Section III. Section IV reports the experimental results and analyses. Finally, Section V concludes this article.

II. METHODS
The two-stage rapeseed mapping process of this work is shown in Fig. 1. In the first stage, rapeseed is identified using the data of the flowering period by CT and VIs constraints. Especially, when the image of the flowering period is at the peak flowering period, VIs are used to obtain the rapeseed distribution directly. On the contrary, CT is performed on the image to be classified and the image of the peak flowering period from 2020 at first. In the second stage, various classifiers are employed to optimize the initial detection result, including three types of deep learning methods, two machine learning classifiers, and NDYI-based mapping methods. Based on the results of the two-stage rapeseed mapping framework, the performance of the aforementioned classifiers has been demonstrated by comparing the mapping precision. Moreover, the effectiveness of CT is researched and the influence of time-series information is also analyzed in this work.

A. Color Transfer and VI-Based Rapeseed Detection
Acquiring substantial training samples is laborious and even impractical in many applications. To overcome this difficulty, CT-and VIs-based rapeseed unsupervised detection method is proposed.
As mentioned before, NDYI expressed by (1) is effective for representing the yellowness of rapeseed, it can be used to detect the rapeseed distribution in the flowering period. Therefore, the NDYI is calculated and combined with a hard threshold α to detect the rapeseed pixels with high confidence. The value of α is decided by experience. Besides, another VI, NDVI, expressed by (2) is also utilized to enrich the initial detection result. In sum, the specific rules for rapeseed detection are as follows.
In this stage, CT is an optional process due to the difference in the images of the flowering period. The whole flowering stage lasted about 30 days. During this time, the available images of different years are not always at the peak flowering stage. In this situation, the yellowness of the flower is not apparent and varies obviously in different regions. The number of samples extracted by VIs will be adversely affected. Actually, in the four-year flowering data, only the image of 2020 is at the peak flowering period. To enhance the yellowness of rapeseed that are not at the peak flowering stage, the color transfer is utilized to alleviate the reflectance differences of the images at the flowering stage. Concretely, a more general color correction method named CTBI [21] is selected to change the color of images, because it is simple and effective. CTBI [21] borrows one image's color features from another image and mainly includes the following steps.
Suppose that the image to be processed is the source image, and the image with the color style we want to convert is the target image. First, the RGB image is projected into color space lαβ through (3) and (4), which is a reversible process. Althrough this decorrelated color space projection, three color channels can be processed separately.
Second, the means and standard deviations of lαβ are computed along each channel. Finally, lαβ derived from the source image are transformed according to the standard deviations of the target image. Specifically, pixels are subtracted by the mean value of the corresponding channel and scaled with factors determined by the standard deviations. Then, the average of the target image is added to each pixel. Take l channel of the source image as an example, the whole transformation can be expressed as where l s , l t , σ l s , and σ l t represent the mean value and the standard deviation of l derived from the source image and target image, respectively. After processing lαβ channel separately, the result is converted back to RGB space.
The motivation of using CTBI is that CTBI is a simple and effective color changing method. The motivation for using CTBI is that it is a simple and effective color changing method. More importantly, it has no effect on the object's properties other than the color. Moreover, the study regions of different years are the same, most land covers are consistent. Changing the color of the images where the flowers are not fully bloomed can make them look more consistent with the image of the peak flowering stage.

B. Result Optimization With Supervised Classifiers
To fully detect the distribution of rapeseed, five commonly used classifiers are adopted to optimize the rapeseed initial detection result. Concretely, the detection result can provide training and validation samples, which means the initial result can be further optimized by supervised classifiers. The studied classification methods including three types of deep learning methods, 1-D CNN (1DCNN), 2-D CNN (2DCNN), and 3-D CNN (3DCNN) are considered in this work. SVM and RF, two widely used machine learning classification methods are also employed. Moreover, NDYI [calculated as (2)] with a hard threshold is used to identify rapeseed for comparison, which is denoted as NDYI thre .
SVM and RF are often employed as baseline machine learning methods in classification tasks because of their robustness [22]. SVM is a supervised classifier that can distinguish different objects through mapping the input vectors to a high-dimensional feature space and constructing the optimal separating hyperplanes [23], [24]. RF is a combination of the decision tree, and the prediction result is obtained by the majority vote on the output categories of each tree [25], [26]. These two methods have been extensively employed in remote sensing applications such as crop mapping [27], [28], [29].
The methods based on deep learning have strong feature extraction ability, which leads to a series of deep neural networks being proposed for classification, especially CNN [30], [31]. Generally, the CNN consists of multiple convolution filters with learnable weights and biases, which can extract hierarchical contextual features of an image. The major operations of the CNN can be explained as where W t is a convolution kernel with the learned weight of the tth layer, and X t−1 denotes the input feature map. This process includes using convolution operation * to W t and X t−1 with the addition of the bias b t , then adopting an activation function σ outside the convolutional layer. Finally, pooling operation pool p with a p × p window size is performed to obtain the next feature map X t . In this work, 1DCNN, 2DCNN, and 3DCNN classifiers with similar multilayer network architectures are constructed to map rapeseed. Three CNN-based networks are composed of three convolutional layers, three batch normalization layers, three activation functions, one max-pooling layer, and one fully connected layer. For 1DCNN, pixel-wise time sequence vectors are extracted as the input. For 2DCNN and 3DCNN, image patches with the target pixel as the center are extracted as the input. Considering rapeseed detection is a binary classification problem, the sigmoid function is selected as the activation function of the output layer. Three types of CNN models are trained with the Adam algorithm and executed in a PyTorch and CUDA environment. The training process is conducted over 300 epochs. Besides, the initial learning rate is set to 0.0001 and the batch size is set to 256. Binary cross-entropy loss with good generality is adopted as the loss function, which is suitable for binary classification problems.
The input time-series MSIs can be regarded as a tensor I ∈ R M×N×B , M, N, B indicate the width, height, and band number of the MSIs, respectively. Data I was clipped into vectors of size 1 × b and patches of size m × m × b for the training and testing, where b represents the number of spectral channels of the MSIs, equal to B here. m × m is the length and width of image patch. In this work, m was set to 5. Moreover, the three convolutional layers of the CNN-based models learn 64, 128, and 256 convolution kernels. The convolution kernel size was consistent with the spatial and spectral axis. Especially, kernel size of 1DCNN, 2DCNN, and 3DCNN was 1×3, 3×3, and 3×3×3.

A. Study Area
As shown in Fig. 2, the study area is located in the East Dongting Lake, which is the northern region of Hunan Province, China. The longitude range of this area is between 112 • 10 and 112 • 60 E, latitude range is between 29 • 0 and 29 • 30 N. The study area mainly contains an alluvial plain with an elevation lower than 50 m. The area covers a total area of 3210 km 2 and mainly includes the north part of Yiyang city and the west part of Yueyang. Influenced by the subtropical monsoon climate, the Dongting Lake area is prone to rainy and cloudy weather in winter and spring. Thus, the available optical remote sensing data are limited. As one of the commodity grain bases in China, the mainly planted crops are paddy rice and rapeseed. Moreover, Fig. 3. Phenological periods of rapeseed in the study area. F, M, and L represent the first, middle, and last ten-day period of a month, respectively.  Fig. 3 presents the phenological periods of rapeseed in the study area.

B. Datasets
As mentioned before, a considerable portion of remote sensing images in the study area are disturbed by the cloudy and rainy conditions. To get enough available images, MSIs captured by multisensor with high quality are collected in this work to explore the robustness of different methods. In particular, five kinds of optical satellites (GF-1, S2, L7, L8, and GF-6) data are utilized. Chinese satellite GF-1 and GF-6 WFV data are downloaded from China Centre for Resources Satellite Data and Application (CRESDA), 1 which can provide red, green, blue (RGB) and near-infrared (NIR) bands with a 16-m spatial resolution and large swath width. Besides, the spatial resolution of RGB and NIR images derived from S2, Landsat satellites are 10 and 30 m, respectively. Table I also lists the detailed information of experimental data. Although the image at the flowering stage is seriously disturbed by clouds in 2019, it is still retained since there is no alternative high-quality MSI can be found.
It should be noted that to keep the temporal information similar in different years, four images of different phenological periods are collected for each year. Besides, there are some differences in band numbers among WFV, Multi Spectral Instrument, ETM+, and OLI sensors. To construct a comparable time series of different years, only the same four bands (i.e., RGB and NIR) of each image are employed in this work. Moreover, Turkoglu et al. [32] also choose the four bands of S2 to realize crop mapping. Because they found that adding more channels increases computing costs, but does not significantly improve performance. 1 [Online]. Available: http://www.cresda.com/EN/ The preprocessing of original GF-1 and GF-6 WFV Level 1 A data includes radiometric calibration, atmospheric correction, and geometric correction. The three correction steps are achieved by the Environment for Visualizing Images (ENVI) 5.3 software. S2 Level 1 C data are preprocessed by Sen2Cor 2 to generate the Level 2 A product. Then, image with RGB and NIR bands is resampled to a 16-m spatial resolution. L7 and L8 data are Collection 2 Level 2 data derived from the United States Geological Survey (USGS), 3 basic preprocessing has been finished in Collection 2. But due to scan-line corrector-off gaps, gap filling is performed on L7 data with ENVI. The processed Landsat data are also resampled to 16-m spatial resolution. Furthermore, all images are registered based on the S2 image obtained on April 12, 2018, because this image have higher resolution and better quality. After registration, all processed images are cropped with a spatial size of 3300×3800 for the subsequent experiment.

C. Training and Testing Samples
In this study, rapeseed mapping is an unsupervised process. In the second stage of our framework, the training samples are automatically generated for supervised classifiers, which includes two categories: rapeseed and nonrapeseed. Moreover, test samples are collected to evaluated the performance of different methods. Since the yellow flower can highlight rapeseed, which is helpful to distinguish rapeseed through human interpretation. Test samples are mainly collected by visual interpretation, and supplemented by ground surveys and Google Earth images. The pixel number of rapeseed test samples from different years   Fig. 4 presents the global and local ground truth in 2017.

D. Validation of Rapeseed Mapping Accuracy
To verify the accuracy of different rapeseed mapping methods, the overall accuracy (OA) and Kappa coefficient are calculated as indicators in this study. OA is the ratio of correctly mapped pixels to all verified pixels, which is suitable for evaluating the overall performance. Kappa coefficient is computed to indicate the consistency of mapping results. Additionally, Precision, Recall, and F1-score based on the confusion matrix (as constructed in Table II) are also used to represent mapping capability. Generally speaking, Precision and Recall are contrary, when one of them is higher, the other often has a lower value. F1-score is the harmonic mean of Recall and Precision, values range from 0 to 1. Considering this article studies a binary classification problem, the specific calculation formulas of the aforementioned indicators are described as follows: where TP, TN, FP, and FN represent the number of true positives, true negatives, false positives, and false negatives cases in the result, respectively. p o is the relative observed agreement among raters (same as OA here), and p e is the hypothetical probability of chance agreement [5]. Besides, the aforementioned evaluation indicators with larger values denote better results.

IV. RESULTS AND ANALYSIS
In this part, the performance of NDYI thre , RF, SVM, 1DCNN, 2DCNN, and 3DCNN are analyzed. The threshold of NDYI thre in this work is decided based on the rule that extracts rapeseed as much as possible without obvious visual misclassification because there is no standard method to find an adaptive threshold. Concretely, the threshold is 0.095, 0.09, 0.08, and 0.09 for 2017-2020 years. Besides, the tree number of RF is 500. The parameters of the SVM are obtained by cross validation with Libsvm. The structure and parameters of CNN models have been introduced in Section II-B. Moreover, parameters without description are set to the default values.
The influence of CTBI on the detection result is also investigated in subsequent experiments. Actually, the rapeseed pixels are determined by NDYI and NDVI at the flowering stage. If the flowers are not in full bloom, CTBI will be performed on the image at first. In addition, threshold α of NDYI is a relatively high value to ensure the reliability of the extracted rapeseed. More in detail, the value of α is 0.13, 0.13, 0.13, and 0.15 for 2017-2020 years, respectively.
Since the initial rapeseed detection result provides a lot of samples, only 10 000 rapeseed samples and 10 000 negative samples are randomly selected for training. The same number of samples is also used for the CNN model validation. When the model gets the best accuracy on validation data, the model is successively adopted to identify rapeseed.
Moreover, to explore the most effective data configuration for rapeseed mapping, the effects of temporal information are analyzed. Multispectral images of different periods are arranged and combined as the input.

A. Effectiveness of Color Transfer
As mentioned in Section II-A, the initial rapeseed detection result of the 2017-2019 years are obtained by NDYI and NDVI after CTBI [21]. In this experiment, multitemporal data from 2017 are utilized to verify the effectiveness of CTBI. Fig. 5 presents the original images and corresponding images after CTBI. As depicted in Fig. 5, CTBI can make the color style of the images consistent with the image from 2020. It should be noted that Fig. 5 shows the normalized original images, which are the processing objects. In other figures, RGB components of MSIs are linearly stretched to make the ground visible. Besides, Fig. 6 displays the NDYI images and corresponding histograms of the original MSIs and the MSIs after CTBI. It can be seen that the distribution of histograms is consistent each year, which means that CTBI does not change data distribution but highlight the rapeseed. In other words, CTBI can help to distinguish rapeseed from other objects.
To understand the effect of CTBI further, a machine learning method (SVM) and a deep learning method (2DCNN) are utilized to identify rapeseed. The number of positive and negative training samples is the same, with a total of 20 000. In addition, the threshold α is 0.11 to get more samples when not performing CTBI. Fig. 7 and Table III present the rapeseed mapping and quantitative results with CTBI or not. According to Fig. 7, a large area is mistakenly divided into rapeseed without CTBI.    Therefore, the value of Recall is higher. In other indicators, the result with CTBI is much higher. For example, OA increased from 0.4960 to 0.9657 and F1-score increased from 0.4212 to 0.8995 without CTBI in SVM results. In terms of the visual performance or evaluation indexes, the results with CTBI are significantly superior. Through this experiment, it can be concluded that color transfer can improve the mapping performance when flowers are not in full bloom.

B. Mapping Accuracy of Different Classifiers
Multitemporal MSIs (16 bands total per year) are used to train and test six different mapping methods mentioned before for the study area. The mapping accuracies that are calculated as the average values of the ten results with random samples (except NDYI thre ) are presented in Table IV. The corresponding rapeseed mapping results are shown in Fig. 8. In addition, four different areas in the red boxes in the Fig. 8 are selected to display the local results, as presented in Fig. 9.
From Table IV, it can be observed that NDYI with threshold has obtained the highest Precision value for the four years, which means the detected rapeseed pixels are more reliable. However, NDYI thre cannot detect all rapeseed pixels with a hard threshold due to phenological differences. The performance of machine learning methods SVM is better than RF. The SVM also can obtain higher evaluation indexes than the deep learning method 1DCNN, especially in 2018 and 2019. For CNN classifiers, 3DCNN gets the highest values of Recall, Kappa, OA, and F1-score in 2017 and 2020. 2DCNN presents the highest values of Kappa, OA, and F1-score in 2018 and 2019. The patch-based 2DCNN and 3DCNN achieve a better mapping accuracy than pixel-based 1DCNN and other methods, since more spatial information is considered.
According to Fig. 8, influenced by the cloud, some regions with high radiation values are misclassified in 2019. Moreover, the wrong classification areas on the rapeseed map are mainly located in the wetlands near the Dongting Lake. One reason is the characteristics of wetland plants are similar to that of rapeseed. Wetland plants are green in spring and sometimes submerged by rising water when rapeseed is harvested. Besides, another difference between the results of different methods is the integrity of the farmland with rapeseed. In order to see the difference more clearly, local areas of the exhibition include the wetland of Dongting Lake and the rapeseed concentrated planting area.
From the first row of Fig. 9, it can be seen that the green plants are prone to misclassification. Objects with high radiation values are also easy to confuse, such as the road and land in the second row of Fig. 9. In the rapeseed planting concentrated area (the third and fourth rows of Fig. 9), the results of 2DCNN have a good spatial agreement with the images from a visual interpretation. Moreover, the interior of the fields is more accurate. According to the aforementioned experiments, 2DCNN achieved better results since it can obtain relatively higher accuracy and visual interpretation in terms of local and global performance. Moreover, the gap filling of Landsat 7 has effects on the result. The pixel value of the filled area is not necessarily true, which also results in some misclassification. When higher quality data are available, this phenomenon will not exist.

C. Performance Analysis of Temporal Information
In order to compare the effects of time-series information on the rapeseed mapping task, images in different periods are combined and used as the input of two classifiers (SVM and 2DCNN). Concretely, four different time feature configurations are designed as Table V shown. Fig. 10 gives the rapeseed maps of all time sequences. Besides, 20 000 training samples are randomly selected to evaluate the effects of different time sequences in this experiment.
According to Table V and Fig. 10, the following conclusions can be drawn. First, the image captured during the flowering stage is crucial for identifying rapeseed. Based on the result of Sequence 1 and Sequence 2 with the SVM classifier, using one image as input can even obtain a better result than using all available MSIs without the flowering period, which illustrates the importance of flowering data for identifying rapeseed. For 2DCNN, Sequence 2 is able to obtain a relatively higher accuracy than Sequence 1. However, when adding data of flowering period and reducing data of repeated phenological periods, OA of rapeseed can be improved from 0.9622 (Sequence 2) to 0.9684 (Sequence 3). Furthermore, as the blue boxes presented in Fig. 10, the result map of Sequence 2 exists obvious misclassification regions. Some vegetation is mistaken for rapeseed. According to the map of Sequences 3 and 4, the results are more accurate when utilizing the information from the flowering period.
Second, four MSIs of different periods are enough for rapeseed mapping. Comparing the results of using four MSIs (Sequence 3) as input with all available MSIs (Sequence 4) as input, the latter is slightly better. In particular, 2DCNN with four MSIs can achieve the best Precision. OA and F1-score of Sequence 4 are only 0.0001 higher than Sequence 3 with the SVM classifier. OA of Sequence 4 is 0.006 higher than Sequence 3 with the 2DCNN classifier. The value of other evaluation indicators and visual performance is also similar, which represents that adding more temporal information about the same phenological periods does not improve the accuracy significantly, but increases the computation cost.

D. Computing Cost
In this article, SVM and RF are realized on the 3.7-GHz CPU and 32-GB memory computer with MATLAB 2018b. CNN models are executed in a PyTorch with CUDA environment, using a GTX 1080 Ti GPU.      testing time of the aforementioned studied classifiers for 2017, where the data preprocessing time is not considered. As Table VI shown, NDYI thre takes the shortest processing time since it does not need learning. Besides, 3DCNN with more parameters needs the longest time for training and testing. Relatively speaking, the time of all methods is acceptable except 3DCNN.

V. CONCLUSION
This work proposes a framework for rapeseed mapping, which utilize the phenological characteristic of rapeseed to realize unsupervised mapping. Specifically, NDYI and NDVI calculated by the flowering period multispectral image are adopted to get an initial rapeseed distribution map. If the flowers are not fully bloomed, a color transformer method CTBI is performed at first. Besides, several universally classifiers are utilized to enrich rapeseed map and compared in this article. More importantly, we analyze the rapeseed mapping task from three perspectives, including the effectiveness of the rapeseed rough detection method with the color transfer, the performance of different classifiers, and the role of temporal information. Experiments conducted on the multitemporal and multisensor MSIs in the Dongting Lake area demonstrate that the framework can identify rapeseed without human labeled samples. Besides, several conclusions can be summarized. 1) Color transformer is able to improve the final mapping accuracy by getting more samples when the rapeseed flowers are not fully bloomed. 2) Considering the mapping accuracy, visual performance, and processing time, 2DCNN is the most appropriate method for rapeseed mapping. 3) Repeated phenological information cannot improve the performance of rapeseed mapping significantly. These findings are reliable since they are obtained based on four years of experiments in the study region covering a total area of more than 3000 km 2 . In other regions similar to the study region, these findings are also reliable in theory.
Moreover, the rapeseed mapping performance can be further improved by optimizing the initial training samples or developing a more effective classifier which can extract discriminant features. This method also can be utilized in more regions with enough available remote sensing data.