Application of Bispectral Analysis to Assess the Effect of Drought on the Photosynthetic Activity of Lettuce Plants Lactuca sativa L.

Abstract

This article proposes a new method for determining the pathological state of a plant, based on a combination of the method for measuring the dynamics of photosystem II pigment fluorescence in the leaves of L. sativa plants and analyzing the resulting time series using bispectral analysis based on the wavelet transform. The article theoretically shows a possible mechanism for the appearance of a peak on the map of bispectrum indexes during nonlinear analog conversion of a physiological signal in a biological object. The phenomenon of increasing the degree of nonlinearity in the transmission of an external periodic signal in plant signaling systems has been experimentally demonstrated.

1. Introduction

The state of a plant’s photosynthetic apparatus affects the plant’s productivity due to its impact on the efficiency of converting energy from visible photons into chemical bond energy during the synthesis of carbohydrates. Among the non-invasive continuous methods for recording the physiological state of a plant, one of the most informative is the method of pulsed fluorimetry of photosynthetic pigments [1,2,3]. This method is capable of determining the state of the plant’s photosynthetic apparatus, the degree of illumination of the plant, and the reserve capabilities of the photosynthetic system.

The theory of amplitude analysis of fluorescence change curves of photosystem pigments includes: (1) measurement of the base level of fluorescence in a plant fully adapted to darkness (F0); (2) measurement of the peak fluorescence value when a dark-adapted photosystem is exposed to a short-term powerful saturating flash of light, the spectral characteristics of which correspond to the absorption spectrum of pigments (Fm); (3) measurement of the peak fluorescence value when a light-adapted photosystem is exposed to a short-term powerful saturating flash of light, the spectral characteristics of which correspond to the absorption spectrum of pigments (F′m); and (4) measuring the value of the basic level of fluorescence in a plant immediately after a saturating flash (F′0) [4]. The algebraic relationships of these parameters make it possible to calculate various parameters describing the physiological state of the plant (maximum quantum yield of photosystem-mediated electron transfer, potential activity of the photosystem, non-photochemical quenching, and linear and cyclic electron flows) [1,4,5]. The use of a PAM fluorimeter makes it possible to carry out both the measurement of the above fluorescent parameters in an automated mode, and the analysis of fluorescent signal records in order to calculate parameters describing the physiological state of the plant. It is precisely this methodological and instrumental approach that is demonstrated by the overwhelming number of works using pulsed fluorimetry of photosystem pigments [6,7]. Dynamic studies of these parameters are carried out quite rarely; in all the works we found, there were typical limitations associated with the inability to establish which features led to certain results. For example, in [8], the wavelet spectra of the initial recordings of a PAM fluorimeter with a duration of 800 s were studied. The recordings included the processes of adaptation of the plant leaf to the light regime (saturating flashes of light); therefore, it is difficult to establish which features of the fluorescent responses led to the results described in the article. Of the listed parameters of PAM fluorimetry, the basic level of fluorescence in the dark mode (F0) is suitable for continuous monitoring of long-term changes, because its measurement can be carried out with a sufficiently high sampling rate, and at the same time, this parameter can be recorded for a long time under stable conditions, which is a necessary for the analysis of the slow dynamics of oscillations.

Previously, we studied low-frequency oscillations in various physiological systems [9,10,11,12] using wavelet transforms. Specifically, we used continuous complex Morlet wavelet transformation [10]. However, recently, many new options have appeared; for example, wavelet with log-normal frequency response. This option has several advantages, both in terms of the quality of the conversion result and the required computational costs. The calculation of a log-normal wavelet occurs by converting the signal into a frequency representation using a fast Fourier transform (FFT), multiplying the Fourier transform by the frequency response of the wavelet mother function, and then applying an inverse FFT to obtain the wavelet coefficients. Since optimized libraries have been developed for FFT using parallel computing, this approach turns out to be less computationally expensive. The log-normal frequency response of a given wavelet transform filter also has the advantage of being symmetrical on the logarithmic scale of frequency or wavelength, i.e., the same shape of the amplitude characteristic both in frequency representation and when displayed in units of wavelengths. In addition, this transformation (through frequency mapping) for all frequencies under study produces a series of wavelet coefficients of the same length without vignetting along the boundaries, i.e., it automatically gives “adaptive” properties [13] without additional algorithmic tricks.

At the moment, there are several analytical approaches to second- and third-order spectral analysis based on Fourier and wavelet decomposition, for example, phase coherence analysis [14] and bispectral analysis [15]. These methods make it possible to analyze the relationship between complex frequency components of signals, and thus take into account the phase difference, both when analyzing synchronized two-channel recordings at the same frequencies (in the case of phase coherence), and cross-frequency relationships in both single-channel and two- and even three-channel recordings (in the case of bispectral analysis). Since we have one channel available to us when recording a fluorescent response, when analyzing inter-frequency interactions, we also consider the results of the bispectral analysis of fluorescence recordings.

However, there are only a relatively small number of studies that have studied the long-term dynamics of fluorescence parameters [16], and even fewer studies using adequate methods of time-frequency spectral analysis. For example, in the work of Shurygin et al. [17], the authors carried out a time-frequency analysis of chlorophyll fluorescence parameters in apple trees with a homemade pulsed fluorometer using high pulse repetition rate technology and the Morlet wavelet transform. They studied the dynamics of chlorophyll fluorescence in plant branches over records about 14 months in length. The authors obtained the maximum wavelet coefficients of the quantum yield of photosystem II at periods of 24 h and 12 h, which indicates the nonlinear nature of the response of this parameter to setting the daily rhythm at 24 h. In another work [8], the authors implemented a fairly straightforward approach to the analysis of fluorescence records for chlorophyll obtained during the standard protocol for determining a set of parameters (Fm, F′m, F0, F′0). Since these recordings represented a response to numerous impulse influences under conditions of plant adaptation to growth lighting, i.e., a highly nonstationary dynamic process, the authors obtained a huge range of wavelet coefficient values using wavelet maps. A further attempt to analyze variations in wavelet coefficients (essentially, a convolution of the moduli of wavelet coefficients, if the authors had used a complex transformation) allowed the authors to distinguish changes in the state of the plant at different times during a sunny day. In [18], the authors studied the parameters of chlorophyll fluorescence in potato leaves based on continuous wavelet transform and spectral analysis; however, they used a mapping pulse fluorimeter, and it was static two-dimensional images of leaves that were subjected to two-dimensional wavelet transform.

The development of methods for analyzing data from continuous monitoring of the physiological state of a plant remains an important scientific task. In the future, this will make it possible to solve the problem of developing a fully automated greenhouse management system for growing vegetable and fruit crops in greenhouses, by tying feedback from the system for assessing the physiological state of the plant to actuators (lighting, watering, nutrients) in a self-regulating system for maintaining optimal growing conditions for plants.

Thus, the aim of this work is to develop a method for studying the long-term dynamics of physiologically important plant parameters: the level of fluorescence of pigments of photosystem II, with periodic test exposure in order to identify possible changes in signal characteristics under drought conditions as a stress factor. The possibility of using bispectral analysis as an applied tool for assessing the physiological state of plants by processing long-term records of the fluorescence of pigments of photosystem II is investigated.

2. Materials and Methods

2.1. Plant Object

Lettuce plants (Lactúca satíva L.) were used as samples. Lettuce sprouts were grown in pots (V = 0.3 L). Peat soil for seedlings was used as a substrate. The height of the plants used for the experiments was 25 cm. Before the experiments, the plants were grown under a 16 h/8 h (day/night) light regime and a temperature of 22–23 °C during the day and 16–17 °C at night, with a humidity of 65%. Lighting was provided by a combination of LEDs (JH-5WBVG14G24-Y6C, Ledguhon, Guangzhou Juhong Optoelectronics Co., Ltd., Guangzhou, China), incandescent lamps (40 W, Electrolighting LLC, Tver, Russia), and UV fluorescent lamps (Litarc Lighting and Electronic Ltd., Shenzhen, China). Figure 1 shows representative photographs of lettuce plants before the drought modeling experiment (Figure 1a) and after (Figure 1b).

2.2. Experimental Design

The plants were placed in an experimental chamber without access to external light and adapted for 24 h. The growth light was not turned on during the entire experiment to exclude the imposition of a circadian rhythm. To measure chlorophyll fluorescence, DUAL-PAM-100 (Waltz, Eichenring, Effeltrich, Germany) was used. The period of switching on the saturating light (300 ms pulse duration, λ = 625 nm, 12,000 µmol photons s⁻¹ m⁻²) was 1 h, and the fluorescence measurement period was 11 s. Measurements were carried out on whole leaves, without separating them from the plant, at room temperature (25 °C) and a humidity of 65%. The experiments were carried out in a darkened room. The experiment with one sample was carried out for 6 days in two stages. During the first 3 days, the plant was provided with the necessary watering. At the same time, the base level of fluorescence was recorded, with a measurement period of 11 s. During the next 3 days, the plant was not provided with water, and the base level of fluorescence was also recorded with a measurement period of 11 s. Figure 2 shows a visualization of the experiment. As a result, two records comprising about 23,000 measurements were obtained from each plant.

2.3. Bispectral Analysis Method

The bispectral analysis method based on the wavelet transform, or wavelet bispectrum [15], is related to second-order spectral methods, similar to the wavelet coherence method [14]; that is, it involves a two-stage mathematical transformation. At the first stage, the signals under study are transformed from a time representation into a frequency-time complex representation, using a complex wavelet transform:

$$\Phi (t,\omega )={\displaystyle \underset{-\infty }{\overset{\infty }{\int }}f(t_i)\Psi (t_i-t,\omega )d{t}_i,}$$

(1)

where f(t) is the original signal or dynamic series, Ψ(t_i − t, ω) is the wavelet filter in the time representation, and Φ(t, ω) is the wavelet coefficient. In this way, information is obtained from the signals on how the amplitudes of the rhythmic components and their phases at different frequencies have changed over time. That is, a two-dimensional map is obtained for the signal, on which time is represented along one axis (usually the abscissa axis), and frequency or period along the other axis (usually the ordinate axis), and at each point of the map the parameter of the oscillatory process is displayed: for example, the modulus of the wavelet coefficient.

At the second stage, the components of the frequency-time representation (complex wavelet coefficients) are subjected to an integral transformation in order to identify correlations in the amplitudes and phases of the rhythmic components. The construction of the integral transformation determines which correlations are identified. In the case of bispectral analysis of a single signal, this integral transformation has the form:

$$\Phi (t,\omega )=A(t,\omega )e^{j\phi },\phi =(\omega t+{\varphi }_1),$$

(2)

$$\rho ({\omega }_1,{\omega }_2)={\displaystyle \underset{-\infty }{\overset{\infty }{\int }}\Phi (t,{\omega }_1)\Phi (t,{\omega }_2)}\overline{\Phi (t,{\omega }_1+{\omega }_2)}dt=$$

(3)

$$={\displaystyle \underset{-\infty }{\overset{\infty }{\int }}A(t,{\omega }_1)}A(t,{\omega }_2)A(t,{\omega }_1+{\omega }_2)e^{j({\phi }_1+{\phi }_2-{\phi }_3)}dt,$$

(4)

where Φ(t, ω) is a complex wavelet coefficient depending on the moment of time and frequency, A is its modulus, ϕ is its phase, and ρ is a bispectral index depending on the frequencies ω₁ and ω₂, which is a time convolution of the product of wavelet coefficients for the frequencies ω₁, ω₂ and ω₁ + ω₂. Such a formulation of the bispectral index simultaneously solves two problems: firstly, it allows us to zero out the wave part of the complex phase of the product, and secondly, it allows us to study the correlation of rhythmic components at two frequencies (ω₁, ω₂), as well as at the frequency of the sum ω₁ + ω₂. In other words, to obtain a non-zero bispectrum index ρ(ω₁, ω₂), the following conditions must be met:

• the signal must simultaneously contain non-zero moduli of wavelet coefficients at frequencies ω₁, ω₂ and ω₁ + ω₂,
• the phase shift of these wavelet coefficients must also be constant.

2.4. Modeling of Analog Nonlinear Transformation of a Harmonic Signal

In order to understand where a signal in which non-zero bispectral indexes appear can appear in a natural system, let us consider a simpler case of studying a bispectral index for coinciding frequencies, i.e., ω₁ = ω₂ = ω:

$$\rho (\omega ,\omega )={\displaystyle \underset{-\infty }{\overset{\infty }{\int }}\Phi (t,\omega )\Phi (t,\omega )}\overline{\Phi (t,\omega +\omega )}dt=$$

(5)

$$={\displaystyle \underset{-\infty }{\overset{\infty }{\int }}{A}^2(t,\omega )}A(t,2\omega )e^{j(2{\phi }_{\omega }-{\phi }_{2\omega })}dt,$$

(6)

it is clear that a search is being made for frequencies for which, in addition to the fundamental frequency, there is a second harmonic, which, for example, may appear in the signal as a result of a purely analog transformation of a harmonic signal in a highly nonlinear system. To illustrate this process, we will model the following signal:

$$x(t)=\sin (\omega t),$$

(7)

$$y(x(t))=A{x}^2(t)+Bx(t)+C,$$

(8)

where x(t) is a harmonic signal, and y(t) is a function of nonlinear transformation of the original signal x(t). In fact, y(t) models the expansion in a Taylor series up to a quadratic term of some nonlinear transfer function in a biological system. The parameters of the model are given in

Table 1. The parameters of the model (Equation (8)).

	Nonlinearity Level (A/B)
	0	0.1	0.3	1.0
ω, Hz	2.778 × 10−4	2.778 × 10−4	2.778 × 10−4	2.778 × 10−4
A	0	0.140	0.397	0.743
B	1.411	1.409	1.325	0.743
C	−0.002	-0.072	−0.269	−0.605

.

The level of nonlinearity is set as the ratio A/B. The model parameters are selected so that the average signal value over the entire time period is zero, and the root mean square deviation of the signal value (RMS) is the same for all signals and equals one. To simulate noise, uniformly distributed random noise with an amplitude of 10 was added to the signal obtained for each A/B variant. This yielded a signal whose main properties repeat those observed in the baseline plant fluorescence signal.

The result of the wavelet transform of the model signals is shown in Figure 3. At nonlinearity levels of 0 (Figure 3a), 0.1 (Figure 3b), and 0.3 (Figure 3c), only the band corresponding to a frequency of 0.278 mHz (period of 3600 s) is clearly visible, which corresponds to the fundamental harmonic of the model signal. At a nonlinearity level of 1.0 (Figure 3d), a band with a frequency of 0.556 mHz (period of 1800 s) is also observed, i.e., this band corresponds to the second harmonic.

Figure 4a shows the time-averaged wavelet coefficients, i.e., the amplitude-frequency characteristic of the signals, which shows that the second harmonic is noticeably present at nonlinearity levels of 0.3 and 1.0. This is a completely expected result, which can also be obtained analytically: with an increase in the nonlinearity level, the amplitude of the second harmonic of the signal increases. Some decrease in the first harmonic with an increase in the nonlinearity level is associated with the redistribution of the oscillation energy towards the second harmonic. When constructing the model signals, we set the root-mean-square deviation of the signal value (RMS), which is an estimate of the total energy of the signal, to one.

The result of the bispectral analysis of signals for all calculated levels of nonlinearity is shown in Figure 5. Two-dimensional maps ρ(ω₁, ω₂) are presented, i.e., the dependence of the bispectral index module on two frequencies. Since a single-channel signal was analyzed, the formula for the bispectrum is symmetrical with respect to the permutation of frequencies ω₁ and ω₂. Therefore, the map ρ(ω₁, ω₂) will also be symmetrical with respect to the diagonal ω₁ = ω₂, so only one half of the map is presented. The figure shows that a noticeable spot (i.e., large bispectrum indices) is observed only on the diagonal elements (i.e., at ω₁ = ω₂), and only for a frequency of 0.278 mHz (period of 3600 s) for nonlinearity levels of 0.3 (Figure 5c) and 1.0 (Figure 5d). In the latter case, the bispectrum indices are higher. For nonlinearity levels of 0.0 (Figure 5a) and 0.1 (Figure 5b), the bispectrum indices differ little from zero. A more detailed representation of the diagonal elements of the bispectrum indices is shown in Figure 4b. It shows that the bispectrum indices for nonlinearity levels of 0 and 0.1 are 10 times smaller than for nonlinearity level 1.0. In this case, small peaks are also noticeable for the frequency of 0.139 mHz, which is caused by the interference of the high-amplitude harmonic fundamental with noise. This also causes the appearance of a peak for the nonlinearity level of 0 at a frequency of 0.278 mHz.

As a result of the simulation, the relationship between the values of the nonlinearity level of the signal transmission system and the bispectrum indices obtained as a result of signal analysis is shown, i.e., an increase in the nonlinearity level during analog conversion of the original pure harmonic signal leads to an increase in the bispectrum index, and precisely at the frequency of the original harmonic.

2.5. Data Analysis

Data analysis via the bispectral analysis method was performed using a script based on the work [15] in Matlab. The generation of data series for modeling the nonlinear transformation of a harmonic signal and the plotting of graphs were performed using MS Excel. Testing of statistical hypotheses using ANOVA was performed in SigmaPlot 11.

3. Results

Analysis of Fluorescence Recordings of Plant Pigments

The study of the signal transmission pathway in physiological systems using bispectral analysis can be carried out in two ways: either using some existing rhythm through which the signals that are products of this rhythm are studied, or using a rhythm imposed on the physiological system by some external influence; in this case, the signals that are its products are studied. Plants depend on diurnal and seasonal rhythms. The oscillations may have long periods. In order to test the hypothesis about changes in the nonlinear properties of the signal transmission pathway in plants, it was decided to create an experimental system with an imposed rhythm. Flashes of white light (white LEDs, 12,000 μmol/(m² × s)) saturating the electron transport chain were chosen as such an influence. The duration of the saturating flash was 0.5 s, and the flash repetition period was 1 h. The period was chosen based on the possibility of collecting a signal duration sufficient for analysis in 72 h, ensuring sufficient accuracy of the oscillation phase assessment with a period of 1 h. An eight-hour fragment of the initial fluorescence recording, registered from a plant under normal conditions, is shown in Figure 6. In the recording, the moment of reaction to the saturating flash was omitted, i.e., the first 2 s after the saturating flash were not recorded, so only the very end of the reaction was recorded. It is interesting that after the fluorescence level was restored to the base level, a quenching process was observed, i.e., the fluorescence level decreased below the base level for a period of about 3 min under normal conditions and for a period of about 7 min under drought conditions, which is clearly seen in the insets of Figure 6, with reactions to the flash aligned with the beginning of the pulse. The effect of protective fluorescence quenching [14] was obviously modified. No other oscillatory processes were visible in the initial recordings shown in Figure 6.

The result of the complex continuous wavelet transform of the plant fluorescence record under normal conditions is shown in Figure 7a, and the result of the complex continuous wavelet transform of the plant fluorescence record under drought conditions is shown in Figure 7b. The graphs of the wavelet coefficients (Figure 7a,b) of the daily records clearly show the lines corresponding to the period of saturating light flashes of 1 h, as well as the doubled (second harmonic) and tripled (third harmonic) periods. There are no noticeable differences between these graphs. The presence of harmonics indirectly indicates the nonlinearity of the reaction process to flashes, but since the flashes themselves were rectangular pulses containing a set of harmonics in the spectrum, this in itself does not indicate the nonlinear nature of the system of signal transmission to plants. However, there is a significant difference between the maps of the values of bispectral coefficients (bispectral maps). In the case of the healthy plant data analysis (Figure 7c), the bispectral analysis identified only the high-frequency region (in the region of 1000 s) of the oscillation period. However, on the bispectral graph for a plant under drought stress conditions (Figure 7d), a spot appeared corresponding to a period of 3600 s, i.e., one hour, i.e., the rhythm of saturating flashes appeared. We obtained records from several plants (three control or healthy plants, and six plants under stress conditions). For analysis, the bispectral coefficients for coinciding frequencies were identified.

The graph of the dependence of these coefficients on frequency (Figure 8a) clearly shows that for a period of 3600 s (1 h) for healthy plants (solid lines) a dip was observed, and for plants under stress conditions (dashed lines), increased values were observed. Testing the hypothesis about the difference in these values (Figure 8b) yields a reliable difference, with an error probability of p = 0.017. This allows us to state that bispectral analysis enables us to identify the state of stress in plants, and that under stress the plant modifies the fluorescent signal of pigments in a nonlinear manner, which forms a phased second harmonic of the reaction to a saturating flash. This result allows us to use continuous recording of fluorescence to identify pathological or stressful conditions in plants directly in open or closed ground.

4. Discussion

In this work, we combined the capabilities of two methods for the first time: the study of the dynamics of plant photosystem II pigments using pulse fluorimetry [1] and the analysis of a dynamic data series using bispectral analysis based on wavelets [15] in order to analyze the internal properties of the plant physiological signaling system. The method of bispectral analysis based on wavelets is quite new and has thus far been used to analyze processes in electrical engineering [19,20] and wave processes in the ocean [21,22], and as a method for assessing nonlinearities in signal propagation in electronics [23]. Bispectral analysis has the advantage of allowing one to study cross-frequency relationships both in single-channel signals and in synchronously recorded two- and three-channel records, which distinguishes it favorably from wavelet-based coherence analysis [14].

In biological systems of various levels of organization, oscillating signals capable of propagating through a biological object are widely represented [24,25], and in the process of propagation, the signal is transformed, propagating in a nonlinear transmission medium. Therefore, the ability to study a signal at different sites of transduction is an opportunity to study the properties of the transduction system itself. The basis of signal transduction in biological systems is the implementation of a sensory function that has an input chain with ligand–receptor interaction, a dose–effect relationship, which in the simplest case is a sigmoidal relationship, and in the simplest case is based on the Michaelis–Menten plot [26]. This graph has a relatively linear section in the vicinity of EC₅₀, but strongly nonlinear behavior when passing to saturation. Therefore, the nonlinear behavior of signal transduction systems in biological objects is perhaps not a surprising fact.

In this work we solved two problems: we investigated the possibilities of bispectral analysis to detect nonlinear behavior in a simple model single-channel case, and we investigated the change in the properties of a nonlinear transduction chain of a simple external periodic signal in a plant object. To solve the first problem, we formed a model signal consisting of a simple harmonic sine function, the result of its transformation in a nonlinear system with different levels of nonlinearity and random noise added to it. The model parameters were selected in such a way as to ensure the same power of the periodic signal in all cases. The power of the added noise was also the same in all cases (the amplitude of the uniformly distributed noise was five times greater than the amplitude of the periodic signal). In the spectral analysis of the obtained signals, we obtained a clear increase in the power of the second harmonic from the level of nonlinearity of the harmonic signal transformation (Figure 3 and Figure 4a) and a rise in the high-frequency region, which is natural for the contribution of uniformly distributed noise. As a result of the bispectral analysis of the simulated signal, the relationship between the values of the nonlinearity level of the signal transmission system and the bispectral indices obtained as a result of the signal analysis is obvious (Figure 5 and Figure 4b), i.e., an increase in the nonlinearity level during analog conversion of the initially pure harmonic signal leads to an increase in the bispectral index, and precisely at the frequency of the initial harmonic. An amazing property of the bispectral analysis is the filtering of incoherent noise and oscillations. That is, neither noise nor an incoherent harmonic signal (in the case of the nonlinearity level of 0) result in a significant response on the map of the bispectral indices. The observed almost complete coincidence of the graphs of the dependence of the diagonal elements of the bispectral indices on frequency for nonlinearity levels of 0 and 0.1 and the peaks present for periods of 7200, 3600, and 1800 s (Figure 4a) is explained by the reaction of the method to a large noise amplitude, which determines the sensitivity of the method.

Thus, the method of bispectral analysis using wavelet transform is a good test (litmus paper) for the level of nonlinearity of analog signal conversion in the system under study, but for this, the system must have some initial harmonic signal. In the case of a plant object, such a signal can be an existing periodic process in the form of a daily change in the intensity of illumination, but the required length of the time series in this case will exceed ten days, since the longer the recording duration, the higher the accuracy of the oscillation phase difference assessment, which is critical for a method operating with complex components of oscillatory processes. However, with such durations, it is difficult to ensure the constancy of other influencing factors, such as air temperature and humidity, soil moisture, nutrient concentration, and exposure to infectious pathogens. However, it is possible to reduce the exposure period to one hour and completely exclude rhythmic exposure with a period of 24 h. Thus, the required recording time can be decreased by 24 times and sufficient accuracy and stability can be obtained with a recording length of 3 days. To simplify the design of the experimental system, we used the DUAL-PAM-100 built-in saturating radiation source for periodic exposure to the plant object.

Figure 6 shows fragments of the initial records of photosystem II pigment fluorescence. The most interesting thing in the records is the change in fluorescence intensity after the end of the saturating flash. It is possible to identify a modification to fluorescence quenching, which is expressed as a decrease in the rate of recovery of fluorescence intensity to the base level under stress conditions (Figure 6). It is quite possible that in the reaction to the flash, there is also a slower process of relaxation oscillations of small amplitude, masked by noise, but this is not visually discernible. When examining a conventional wavelet spectrum in the form of a time-frequency map (Figure 7a,b), it is evident that the signal of both the control plant and the plant under drought conditions contains a fundamental rhythm with a period of 3600 s, and its second (1800 s) and third (1200 s) harmonics. Moreover, their amplitudes for the control plant and the plant under drought conditions visually coincide. There is nothing surprising in the presence of higher harmonics during signal transduction in a biological system, since, as stated above, the main transfer functions of ligand–receptor interaction are nonlinear. Bispectral analysis of 3-day fluorescence records (Figure 7c,d) shows their fundamental difference: in the case of drought, a spot appears on the bispectral map, the frequency position of which is determined by the imposed rhythm with a period of 3600 s. Moreover, when analyzing the statistics for diagonal elements (Figure 8a), it is evident that the values of the bispectral indices for control plants at the frequency of the external stimulus differ reliably by two orders of magnitude (Figure 8b). Thus, using bispectral analysis, it was possible to determine the pathological state of the plant using a non-invasive method. An interesting question arises surrounding the reason for the appearance of high bispectral indices with the obvious presence of upper harmonics of the imposed signal for both control and damaged plants. A possible answer is related to the non-conservation of the phase difference for the fundamental period and for the second harmonic in the control plants. This may be related to the coincidence of the second harmonic frequency with some internal physiological rhythm of the plant. Such rhythmic processes have indeed been described in the literature: the process of leaf movement, with periods so short that they were called micronutations (periods from 12 to 30 min), was found in the tendrils of green beans in most, but not all individuals [27]. The frequencies of ordinary nutations are in the range of 50 μHz (periods of about 20–300 min) [28], and transpiration oscillations occur over a period spanning from a few minutes to an hour [28]. Since all of the above rhythms are systemic, they can also influence the signaling of photosystem II, and therefore the fluorescent response.

5. Conclusions

In this work, we have demonstrated the capabilities of bispectral analysis for identifying nonlinear behavior of signal transmission chains in biological systems with single-channel recording of physiological dynamics. We have also demonstrated the possibility of detecting stress in plants during drought using a noninvasive method that combines pulse fluorimetry and bispectral analysis.

The approach we propose may be useful for further practical application in greenhouse farms, where precise detection of the physiological state of a plant is required in order to increase crop yield. It is worth noting that the approach we are considering, although non-invasive, required continuous fixation of the leaf blade in the PAM fluorimeter device, which undoubtedly imposes a number of limitations in large-scale application in agricultural technologies. Consequently, an important unsolved problem remains the development of a contactless system for assessing the fluorescence of plant pigments which can assess whole plants, as well as several plant objects simultaneously.