# Arabidopsis Growth Simulation Using Image Processing Technology

Authors: Junmei Zhang; Qiuhong Ke; Ye Tian
Journal: Journal of Electrical and Computer Engineering
 (2014)

---

## Abstract

This paper aims to provide a method to represent the virtualArabidopsis plant at each growth stage. It includes simulating the shape and providing growth parameters. The shape is described with elliptic Fourier descriptors. First, the plant is segmented from the background with the chromatic coordinates. With the segmentation result, the outer boundary series are obtained by using boundary tracking algorithm. The elliptic Fourier analysis is then carried out to extract the coefficients of the contour. The coefficients require less storage than the original contour points and can be used to simulate the shape of the plant. The growth parameters include total area and the number of leaves of the plant. The total area is obtained with the number of the plant pixels and the image calibration result. The number of leaves is derived by detecting the apex of each leaf. It is achieved by using wavelet transform to identify the local maximum of the distance signal between the contour points and the region centroid. Experiment result shows that this method can record the growth stage of Arabidopsis plant with fewer data and provide a visual platform for plant growth research.


---

## Body

## 1. Introduction

As an important model plant,Arabidopsis has been widely used in plant genetic research [1]. The genome sequencing work was completed at the end of 2000 [2]. To analyse the function of genes, the phenotypes also need to be studied. The phenotypes are the observable characteristics of the plants. They are the inevitable results of the genetic variation and the growth environment [3]. The analysis of plant phenotypes can reveal the impacts of the genes and the environment, which is an efficient method to cultivate new crops.The phenotypes ofArabidopsis include the total area, the number of leaves, and the overall contour. The total area and the number of leaves of a plant reflect its growth situation. The changing of these two parameters can also reveal the growth rate. The overall contour is a significant shape parameter. It can affect the morphology, growth orientation, density, and even size of a plant. With the analysis of the overall contour, we can get a better grasp of the morphological features of the plant.However, the number of contour points is quite large because of the high resolution of image. It will require much storage to directly store the contour pixels. Besides, the plant varies at different growth stages, resulting in different numbers of contour pixels. And it will also add difficulties to store the data for their different lengths.So how to represent the shapes of differentArabidopsis plants? Previous research has shown various ways to describe shape variation [4]. What we have done was using elliptic Fourier descriptors (EFDs) [5] to analyse the data in the frequency domain, instead of the spatial domain. It decomposes a curve into a set of harmonically related ellipses [6–9]. EFDs can be easily normalized to represent shapes unaffected by their scale, location, or orientation. It has been widely used for morphological estimation. Menesatti et al. [10] used the Fourier ellipse analysis to extract features of hazelnuts to estimate the shape. Costa et al. [11] extracted the polar signature of the contour of oranges and calculated the coefficients of Fourier harmonics to describe the size, shape, and orientation of each fruit. Similar work has been carried out by other researchers [12–14].In this paper, EFDs will be used to analyse the contour, and eachArabidopsis plant will be represented with much fewer coefficients, rather than the original contour data. This method not only saves much storage, but also can be used to simulate the shape of the plant. The total area and the number of leaves of each plant will also be calculated. The total area is obtained with the number of the plant pixels and the image calibration result. The number of leaves is derived by detecting the apex of each leaf. It is achieved by using wavelet transform to identify the local maximum of the distance signal between the contour points and the region centroid.


## 2. Image Preprocessing

### 2.1. Image Acquisition

At each pot ofArabidopsis plant, a small square of checkerboard, with the side of 12 mm, was used to do image correction and calibration. With the alternating dark and bright grids, the corners of the checkerboard pattern present as very strong features and can be detected easily [15]. The top view of the plant was taken by a camera (Canon EOS 500D), in which the image resolution was set to 3456 × 2304 pixels. The total number of plants is 100. One of the original images is shown in Figure 1(a).Image preprocessing. (a) Original image of anArabidopsis plant, with soil background. (b) Corner detection result. The corners are marked with red dots. (c) Image correction result.
(a)
(b)
(c)


### 2.2. Image Correction and Calibration

It is inevitable that the image is distorted. It needs to be corrected for further geometrical measurement. The distortion can be expressed as(1)
[
x
1
′
x
2
′
x
3
′
]
=
[
h
11
h
12
h
13
h
21
h
22
h
23
h
31
h
32
h
33
]
[
x
1
x
2
x
3
]
,
where (x
1
,
x
2
,
x
3) is the homogeneous coordinate of a point in real world, (x
1
′
,
x
2
′
,
x
3
′) is the coordinate in the image. The matrix H is the transform function. It can be computed with four pairs of points. With the result, the inverse matrix H
-
1 can be obtained and the distortion can be eliminated by convolving it with the image.In this paper, the corners of the checkerboard were detected to compute the transform function. The checkerboard region was first automatically detected by using (2) to recognize the blue border. With blue border region, the checkerboard region was acquired by filling the holes. The Harris corner detection algorithm was then used to detect the corners [16]. Thus the image could be corrected. Figure 1(b) shows the corner detection result and Figure 1(c) shows the image after correction. Since the side of checkerboard is 12 mm, its real size is 144 mm [2]. The real size of each pixel could also be obtained with (3). Consider
(2)
B
-
R
>
0
,
B
-
G
>
0
,
B
>
150
,
where R
,
G, and B are the red, green, and blue values of a pixel. Consider
(3)
A
d
z
=
144
N
b
,
where N
b is the number of checkerboard region pixels and A
d
z is the real size of each pixel.


## 2.1. Image Acquisition

At each pot ofArabidopsis plant, a small square of checkerboard, with the side of 12 mm, was used to do image correction and calibration. With the alternating dark and bright grids, the corners of the checkerboard pattern present as very strong features and can be detected easily [15]. The top view of the plant was taken by a camera (Canon EOS 500D), in which the image resolution was set to 3456 × 2304 pixels. The total number of plants is 100. One of the original images is shown in Figure 1(a).Image preprocessing. (a) Original image of anArabidopsis plant, with soil background. (b) Corner detection result. The corners are marked with red dots. (c) Image correction result.
(a)
(b)
(c)


## 2.2. Image Correction and Calibration

It is inevitable that the image is distorted. It needs to be corrected for further geometrical measurement. The distortion can be expressed as(1)
[
x
1
′
x
2
′
x
3
′
]
=
[
h
11
h
12
h
13
h
21
h
22
h
23
h
31
h
32
h
33
]
[
x
1
x
2
x
3
]
,
where (x
1
,
x
2
,
x
3) is the homogeneous coordinate of a point in real world, (x
1
′
,
x
2
′
,
x
3
′) is the coordinate in the image. The matrix H is the transform function. It can be computed with four pairs of points. With the result, the inverse matrix H
-
1 can be obtained and the distortion can be eliminated by convolving it with the image.In this paper, the corners of the checkerboard were detected to compute the transform function. The checkerboard region was first automatically detected by using (2) to recognize the blue border. With blue border region, the checkerboard region was acquired by filling the holes. The Harris corner detection algorithm was then used to detect the corners [16]. Thus the image could be corrected. Figure 1(b) shows the corner detection result and Figure 1(c) shows the image after correction. Since the side of checkerboard is 12 mm, its real size is 144 mm [2]. The real size of each pixel could also be obtained with (3). Consider
(2)
B
-
R
>
0
,
B
-
G
>
0
,
B
>
150
,
where R
,
G, and B are the red, green, and blue values of a pixel. Consider
(3)
A
d
z
=
144
N
b
,
where N
b is the number of checkerboard region pixels and A
d
z is the real size of each pixel.


## 3. Image Segmentation

The region of green plant in an image was segmented by using Excess Green (ExG) minus Excess Red (ExR) index [16]. It is given as
(4)
ExG
-
ExR
=
2
g
-
r
-
b
-
(
1.4
r
-
g
)
,
where r
,
g, and b are the red, green, and blue chromatic coordinates. For a 24-bit colour image, they are given as
(5)
r
=
R
R
+
G
+
B
,
g
=
G
R
+
G
+
B
,
b
=
B
R
+
G
+
B
,
where R
,
G, and B are the RGB values of a pixel.The binary images were obtained with fixing threshold 0. Figure2(a) shows the segmentation result of the plant. With the segmentation result, the number of plant pixels could be obtained. Combined with the calibration result in preprocessing step, the real area of the plant could be calculated. It is given as
(6)
A
p
=
A
d
z
N
,
where N is the number of pixels of the plant region.Image segmentation. (a) Binary image of segmentation result. (b) Image after filling the holes. (c) Outer contour of the image.
(a)
(b)
(c)The holes inside the region were filled, as shown in Figure2(b). Thus, the outer overall contour, shown in Figure 2(c), could be obtained.


## 4. Contour Reconstruction

### 4.1. Boundary Tracking

After obtaining the plant region, the boundary pixels need to be stored clockwise or counter clockwise. This could be done by using boundary point tracking algorithm [17]. There are three steps to conduct the tracking.Step 1.
Select the pixel of the minimum row coordinate from the pixels of the minimum column value of the region. Let this pixelb
0 as the starting tracking pixel of the boundary. Another pixel c
0 is the west neighbour points of b
0. Store b
0
,
c
0. Let b
=
b
0
,
c
=
c
0.Step 2.
Search the3
×
3 neighbourhood of b in a clockwise direction from c to find the next boundary pixel b
1. In the 8 neighbour point series of b, the background pixel in front of b
1 is recorded as c
1. Store b
1
,
c
1. Let b
=
b
1
,
c
=
c
1.Step 3.
Repeat Step2 until the pixel b is equal to b
0 and the next boundary pixel is b
1.With these steps, the boundary point series is obtained.


### 4.2. Elliptic Fourier Descriptors

In the boundary series, each pixel has two values. The image can be viewed as a complex space and the two coordinates can be regarded as the real part and imaginary part. Therefore, the boundary series can be expressed as a discrete function shown in(7)
c
(
i
τ
)
=
x
(
i
τ
)
+
j
y
(
i
τ
)
i
=
1,2
,
…
,
m
,
x
(
i
τ
) and y
(
i
τ
) define the two coordinate values of the sampling point i. If the boundary is closed, it can be performed by theFourier expansion in the two coordinate directions to calculateFourier coefficients. It is given as
(8)
c
(
i
τ
)
=
a
x
0
2
+
∑
k
=
1
m
/
2
(
a
x
k
cos
⁡
(
k
ω
i
τ
)
+
b
x
k
sin
(
k
ω
i
τ
)
)
+
j
(
a
y
0
2
+
∑
k
=
1
m
/
2
(
a
y
k
cos
⁡
(
k
ω
i
τ
)
+
b
y
k
sin
(
k
ω
i
τ
)
)
)
,
where
(9)
a
x
k
=
2
m
∑
i
=
1
m
x
i
cos
⁡
(
k
ω
i
τ
)
,
b
x
k
=
2
m
∑
i
=
1
m
x
i
sin
(
k
ω
i
τ
)
,
a
y
k
=
2
m
∑
i
=
1
m
y
i
cos
⁡
(
k
ω
i
τ
)
,
b
y
k
=
2
m
∑
i
=
1
m
y
i
sin
(
k
ω
i
τ
)
.
Thus, eachFourier harmonics has one set of coefficients (a
x
k
,
b
x
k
,
a
y
k
,
b
y
k). m is the number of points. T is the period of the boundary. Generally, it is equal to 2π. τ is the sampling period, which is equal to T
/
m. ω defines the fundamental frequency. It is equal to 2
π
/
T. Thus, k
ω
i
τ is equal to 2
π
k
i
/
m. The harmonics number k defines the level of detail of the characterization [18]. The higher the number is, the more precise descriptions the coefficients provide. However, they also contain much more noise. Generally speaking, the curve can be represented with less than m
/
2 coefficients.Figure3 shows reconstructed contour with different harmonics numbers. It can be seen that the contour is reconstructed well with 30 set of coefficients. As the coefficients increase, there are no more new details of the contour but some noise in the image. In other words, 120 coefficients are enough to store the contour of the plant. The data could be used anytime to reconstruct the contour, providing the same representation as the original contour in Figure 2(c). The original contour contains 1163 pixels. Compared with the original number, the coefficients save a lot of storage.Contour reconstruction. (a)–(f) Contour reconstruction with 10, 20, 30, 40, 50, and 60 harmonics numbers, respectively.
(a)
(b)
(c)
(d)
(e)
(f)


## 4.1. Boundary Tracking

After obtaining the plant region, the boundary pixels need to be stored clockwise or counter clockwise. This could be done by using boundary point tracking algorithm [17]. There are three steps to conduct the tracking.Step 1.
Select the pixel of the minimum row coordinate from the pixels of the minimum column value of the region. Let this pixelb
0 as the starting tracking pixel of the boundary. Another pixel c
0 is the west neighbour points of b
0. Store b
0
,
c
0. Let b
=
b
0
,
c
=
c
0.Step 2.
Search the3
×
3 neighbourhood of b in a clockwise direction from c to find the next boundary pixel b
1. In the 8 neighbour point series of b, the background pixel in front of b
1 is recorded as c
1. Store b
1
,
c
1. Let b
=
b
1
,
c
=
c
1.Step 3.
Repeat Step2 until the pixel b is equal to b
0 and the next boundary pixel is b
1.With these steps, the boundary point series is obtained.


## 4.2. Elliptic Fourier Descriptors

In the boundary series, each pixel has two values. The image can be viewed as a complex space and the two coordinates can be regarded as the real part and imaginary part. Therefore, the boundary series can be expressed as a discrete function shown in(7)
c
(
i
τ
)
=
x
(
i
τ
)
+
j
y
(
i
τ
)
i
=
1,2
,
…
,
m
,
x
(
i
τ
) and y
(
i
τ
) define the two coordinate values of the sampling point i. If the boundary is closed, it can be performed by theFourier expansion in the two coordinate directions to calculateFourier coefficients. It is given as
(8)
c
(
i
τ
)
=
a
x
0
2
+
∑
k
=
1
m
/
2
(
a
x
k
cos
⁡
(
k
ω
i
τ
)
+
b
x
k
sin
(
k
ω
i
τ
)
)
+
j
(
a
y
0
2
+
∑
k
=
1
m
/
2
(
a
y
k
cos
⁡
(
k
ω
i
τ
)
+
b
y
k
sin
(
k
ω
i
τ
)
)
)
,
where
(9)
a
x
k
=
2
m
∑
i
=
1
m
x
i
cos
⁡
(
k
ω
i
τ
)
,
b
x
k
=
2
m
∑
i
=
1
m
x
i
sin
(
k
ω
i
τ
)
,
a
y
k
=
2
m
∑
i
=
1
m
y
i
cos
⁡
(
k
ω
i
τ
)
,
b
y
k
=
2
m
∑
i
=
1
m
y
i
sin
(
k
ω
i
τ
)
.
Thus, eachFourier harmonics has one set of coefficients (a
x
k
,
b
x
k
,
a
y
k
,
b
y
k). m is the number of points. T is the period of the boundary. Generally, it is equal to 2π. τ is the sampling period, which is equal to T
/
m. ω defines the fundamental frequency. It is equal to 2
π
/
T. Thus, k
ω
i
τ is equal to 2
π
k
i
/
m. The harmonics number k defines the level of detail of the characterization [18]. The higher the number is, the more precise descriptions the coefficients provide. However, they also contain much more noise. Generally speaking, the curve can be represented with less than m
/
2 coefficients.Figure3 shows reconstructed contour with different harmonics numbers. It can be seen that the contour is reconstructed well with 30 set of coefficients. As the coefficients increase, there are no more new details of the contour but some noise in the image. In other words, 120 coefficients are enough to store the contour of the plant. The data could be used anytime to reconstruct the contour, providing the same representation as the original contour in Figure 2(c). The original contour contains 1163 pixels. Compared with the original number, the coefficients save a lot of storage.Contour reconstruction. (a)–(f) Contour reconstruction with 10, 20, 30, 40, 50, and 60 harmonics numbers, respectively.
(a)
(b)
(c)
(d)
(e)
(f)


## 5. Leaf Apex Recognition

From the overall contour, it can be seen that the distance between the apex of each leaf and the region centroid is a local maximum. In this case, we firstly calculated the distance between each boundary point and the centroid. Thus we obtain a 1D signal. The next step is to identify the local maximum points of the signal.


### 5.1. 1D Distance Curve

The distance between the contour point and the centroid is calculated as(10)
D
i
=
(
x
i
-
x
o
)
2
+
(
y
i
-
y
o
)
2
,
where (x
i
,
y
i) is the coordinates of the contour point i, (x
o
,
y
o) is the values of the centroid point.


### 5.2. Local Maximum Recognition

The original distance curve signal contains much noise. It is impossible to locate the local maximum only by comparing each contour point with the points at its two sides. In this paper, we used wavelet transform to detect the local maximum point.Wavelet transform is often used to detect the singularity of the signal. When the signal contains much change, the wavelet coefficients will be quite large. While on the contrary, the wavelet coefficients will be zero if the point is a local minimum or maximum. Another merit of wavelet analysis is that it provides a decomposition of the original signal into different scales. So it is very convenient to choose the scale according to different demands to analyse the signal.In the experiment, we first expand the signal into2
n length so that it can be performed by wavelet transform. Then the signal was discomposed into 5 levels with Haar wavelet transform. The wavelet coefficients were used to analyse the local maximum. For a point, if the coefficient of its former one is smaller than 0 while that of its latter one is larger than 0, it will be identified as the zero-crossing point of the wavelet coefficients, as well as the local maximum of the distance signal. Figure 4(a) shows the wavelet coefficients in the last level. The zero coefficients have been marked with red dots. Figure 4(b) shows the local maximum of 1D distance signal. Figure 4(c) shows the leaf apex detection result.Local maximum recognition. (a) Wavelet coefficients in the fifth level. The zero coefficients are marked with red dots. (b) Distance curve between the contour point and the centroid. The local maximum points are labelled with red dots. (c) The apex of each leaf.
(a)
(b)
(c)


## 5.1. 1D Distance Curve

The distance between the contour point and the centroid is calculated as(10)
D
i
=
(
x
i
-
x
o
)
2
+
(
y
i
-
y
o
)
2
,
where (x
i
,
y
i) is the coordinates of the contour point i, (x
o
,
y
o) is the values of the centroid point.


## 5.2. Local Maximum Recognition

The original distance curve signal contains much noise. It is impossible to locate the local maximum only by comparing each contour point with the points at its two sides. In this paper, we used wavelet transform to detect the local maximum point.Wavelet transform is often used to detect the singularity of the signal. When the signal contains much change, the wavelet coefficients will be quite large. While on the contrary, the wavelet coefficients will be zero if the point is a local minimum or maximum. Another merit of wavelet analysis is that it provides a decomposition of the original signal into different scales. So it is very convenient to choose the scale according to different demands to analyse the signal.In the experiment, we first expand the signal into2
n length so that it can be performed by wavelet transform. Then the signal was discomposed into 5 levels with Haar wavelet transform. The wavelet coefficients were used to analyse the local maximum. For a point, if the coefficient of its former one is smaller than 0 while that of its latter one is larger than 0, it will be identified as the zero-crossing point of the wavelet coefficients, as well as the local maximum of the distance signal. Figure 4(a) shows the wavelet coefficients in the last level. The zero coefficients have been marked with red dots. Figure 4(b) shows the local maximum of 1D distance signal. Figure 4(c) shows the leaf apex detection result.Local maximum recognition. (a) Wavelet coefficients in the fifth level. The zero coefficients are marked with red dots. (b) Distance curve between the contour point and the centroid. The local maximum points are labelled with red dots. (c) The apex of each leaf.
(a)
(b)
(c)


## 6.Arabidopsis Growth Simulation

The images ofArabidopsis plants were taken every day during the growth period. The growth simulation is based on the obtained data in 5 days. The experiments were run on a standard PC with 2.1 GHz Dual-Core CPU and 2 GB memory. The processing time is on average 0.5 second per image using MATLAB. Figure 5 shows one example of theArabidopsis plants. Figures 5(a)–5(e) show the overall contour reconstruction result. They were realized with (8). The harmonics number is 30. The apex of each leaf is also marked in each image. With the leaf apex, the number of leaves can be obtained. Figure 5(f) compares total area of the plant in each day.Growth simulation. (a)–(e) Contour reconstruction during 5 days. The leaves apexes are marked with red dots. (f) Increase of the total area of the plant.
(a)
(b)
(c)
(d)
(e)
(f)


## 7. Conclusion

This paper provides a method to represent the virtualArabidopsis plant at each growth stage. TheArabidopsis plant is automatically segmented and the overall contour is obtained. The elliptic Fourier analysis is introduced to extract the coefficients of the contour. With this method, all plant contours could be represented and reconstructed with only 120 coefficients. The data is much less than the original number of contour points, which help to save a lot of storage. In addition, it also makes the data storage much convenient with the same length of coefficients rather than different lengths of contour points. The total area of the plant is measured and the number of leaves is also detected by using wavelet transform to identify the local maximum of the contour-centroid distance signal. Experiment shows that this method can record the growth stage ofArabidopsis plant with fewer data and provide a visual platform for plant growth research.



---

*Source: 751927-2014-08-24.xml*