**Nonstationary Flood Frequency Analysis Using Univariate and Bivariate Time-Varying Models Based on GAMLSS**

#### **Ting Zhang 1, Yixuan Wang 1,\*, Bing Wang 2, Senming Tan 3 and Ping Feng 1**

Received: 16 May 2018; Accepted: 14 June 2018; Published: 21 June 2018

**Abstract:** With the changing environment, a number of researches have revealed that the assumption of stationarity of flood sequences is questionable. In this paper, we established univariate and bivariate models to investigate nonstationary flood frequency with distribution parameters changing over time. Flood peak Q and one-day flood volume W1 of the Wangkuai Reservoir catchment were used as basic data. In the univariate model, the log-normal distribution performed best and tended to describe the nonstationarity in both flood peak and volume sequences reasonably well. In the bivariate model, the optimal log-normal distributions were taken as marginal distributions, and copula functions were addressed to construct the dependence structure of Q and W1. The results showed that the Gumbel-Hougaard copula offered the best joint distribution. The most likely events had an undulating behavior similar to the univariate models, and the combination values of flood peak and volume under the same OR-joint and AND-joint exceedance probability both displayed a decreasing trend. Before 1970, the most likely combination values considering the variation of distribution parameters over time were larger than fixed parameters (stationary), while it became the opposite after 1980. The results highlight the necessity of nonstationary flood frequency analysis.

**Keywords:** nonstationarity; univariate model; GAMLSS; bivariate model; copulas

#### **1. Introduction**

Flood frequency analysis is the premise and foundation of water conservancy project planning and construction. The current flood frequency analysis methods usually assume that the flood series satisfies consistency, that is, that the distribution form or the statistical law of the flood sequence is fixed [1]. However, with climate change and intensification of human activities, especially the construction of large-scale water conservancy and water conservation engineering and the urbanization process, the runoff yield and concentration mechanism, and the temporal and spatial distribution of flooding have been changed [2]. This results in the inconsistency of flood series and the unreliability of the frequency obtained from current frequency analysis methods [3]. Therefore, it is of great significance to study the nonstationary flood frequency analysis methods [4].

Existing nonstationary flood frequency analysis methods in literature include mixture distribution methods, conditional probability distribution methods, and time-varying moment methods. The main idea of the nonstationary flood frequency analysis methods based on mixture distribution is that the individuals of the extreme series are not from the same population [5]. That is, the series formed by different hydrological processes does not follow the same distribution, thus it was assumed to consist of several sub-distributions. The nonstationary flood frequency analysis methods based on conditional probability distribution divide the flood into several periods based on the differences in flood formation mechanism, analyze the occurrence probability of annual maximum value in the different periods, and then obtain the probability density function of the extreme series [6].

Different from the mixture distribution models and conditional probability distribution models, models with time-varying moment consider the change of climate and land surface to have resulted in a change in the physical processes and mechanism of flood generation, such that the parameters of the distribution followed by the flood sequence are functions of time rather than constant. Much attention has been paid to time-varying moment models [7]. Vasiliades et al. [8] applied a time-varying moment model based on Generalized Extreme Value (GEV) distribution to analyze nonstationary frequency, through assuming the parameters of GVE distribution to be the functions of time or other factors and conducting a goodness-of-fit test to the model, thus verifying the significance of the nonstationarity of hydrological sequences.

This study addressed models with time as covariate for nonstationary flood frequency analysis based on GAMLSS (Generalized Additive Models for Location, Scale, and Shape) theory. GAMLSS was first proposed by Rigby and Stasinopoulos [9,10]. This model overcomes the limitations of the GLM and GAM models, greatly expands the range of the distribution types, and provides a variety of ways to produce different distributions, including a series of continuous and discrete distributions with high skewness and/or high peak. In addition, the systematic components provide more plentiful content. For example, it can introduce more complex parametric (linear or non-linear), semi parametric, non-parametric, or random-effect terms to establish models between the distribution parameters (mean, variance, etc.) and the explanatory variables. The GAMLSS model has been widely applied in the military, economics, medicine, and other fields [11–13]. Hydrologists have also done many researches using GAMLSS in recent years. Serinaldi and Kilsby [14] used the GAMLSS model to analyze the monthly rainfall data of 6 stations in England and Wales. They found that the model could better describe the characteristics of rainfall series, and had better performance for fitting the relation between extreme rainfall events, rainfall and atmospheric circulation index, and sea surface temperature. López and Francés [15] proposed two methods based on the GAMLSS model to investigate the nonstationary frequency analysis of the annual maximum flood records of 20 Spanish inland rivers. The results illustrate that the nonstationarity of flood series caused by the effects of climate change and human activities can not to be ignored, and GAMLSS provides a convenient and flexible model framework for considering the influence of climate factors and human activities in the analysis of non-stationary flood frequency.

In flood frequency analysis, univariate probability distribution functions are usually used to estimate the occurrence probability or magnitude of the flood peak or volume in a certain region. However, flood events involve more than one characteristic variable such as flood peak discharge, flood volume, flood water level, etc. In order to estimate the probability of flooding, one needs to know not only the high and extreme values of each variable, but also the likelihood of their occurring simultaneously [16]. The main issue of the univariate models is their difficulties in capturing the underlying joint probability among multiple physical processes, and this will lead to underestimation of the associated occurring probability [17]. For instance, if only the rainfall is considered to estimate the flood risk for a catchment, the resultant estimation would be significantly lower than its true risk, when there is a statistically significant dependence between the rainfall on the catchment and the downstream high water levels. To this end, bivariate models are used to address this issue [18].

Copula functions, for which the marginal distribution of each variable is uniform, were adopted in the bivariate model with time as covariate in this study. They are popular in high-dimensional statistical applications because they allow one to easily model and estimate the distribution of random vectors by estimating marginals and copulae separately. They can describe the linear, non-linear, symmetric, and asymmetric relations between variables, and are simple, flexible, and adaptable in application. Therefore, copulas are effective mathematical tools to construct the multivariate joint distribution and correlation between variables. In recent years, they have been widely used in multivariate hydrological frequency analysis. In drought characteristics analysis, Mirabbasi et al. [19] used a copula function to establish the joint probability distribution between drought duration and drought degree. In rainfall frequency analysis, Zhang and Singh [20] adopted copula functions to construct the bivariate joint distribution between rainfall intensity and depth, rainfall intensity and duration, and rainfall depth and duration, respectively. The results were compared with a Gumbel mixture model and a two-dimensional normal transformation distribution model. Fu G. and Butler D. [21] used the copula method to separate the dependence structure of rainfall variables from their marginal distributions, and analyzed the different impacts of dependence structure and marginal distributions on system performance.

This paper takes Wangkuai Reservoir, which is undergoing substantial change in climate, land use/land cover, and increased number of soil-water conservation projects in Daqing River Basin, to construct both univariate and bivariate time-varying moment models for flood frequency analysis based on GAMLSS theory. The inflow flood peak and flood volume time series (1956–2004) of Wangkuai Reservoir were selected as basic data to discuss the nonstationary univariate and peak-volume bivariate joint flood frequency analysis methods. Flood quantiles and the combined values of flood peak and flood volume under certain exceedance probabilities have been worked out. This study aims to provide new ideas and approaches for nonstationary flood frequency analysis method under a changing environment.

#### **2. Study Region and Data**

Wangkuai Reservoir is located in the upstream of Sha River, Daqinghe Catchment (Figure 1). Its construction started in June 1958 and finished in September 1960. The control area of the reservoir is 3770 km2, and the storage capacity is 13.89 × 108 m3. The currently used design floods were calculated with flood data series under the assumption of stationarity. The watershed receives an average precipitation of 626.4 mm annually, mostly in the summer (70–80%). The annual mean temperature is 7.4 ◦C.

**Figure 1.** The study area: Wangkuai Reservoir in the Daqinghe river basin with the Wangkuai Reservoir watershed in the upper-left corner.

Since 1980, a series of water conservation measures have been carried out in the Wangkuai Reservoir catchment, such as closing land for reforestation and returning farmland to forest. Meanwhile, a number of small hydraulic structures have been constructed. These factors have increased the vegetation coverage rate and significantly changed the land surface, which has affected the flood process in this watershed and thus resulted in nonstationarity of flood series, as revealed by many studies [22–24].

Flooding runoff data have been monitored for a period of 49 years from 1956 to 2004, and collected on hourly basis. The data were provided by Hydrology and Water Resources Survey Bureau of Hebei Province. Maximum flood peak Q and maximum one-day flood volume W1 of each year are selected in this work.

#### **3. Methods**

#### *3.1. Generalized Additive Models in Location, Scale, and Shape (GAMLSS)*

In this study, we adopted a general class of regression models which was called the Generalized Additive Models in Location, Scale, and Shape (GAMLSS) to analyze the nonstationary flood frequency. In GAMLSS models, one assumes that the vector of actual observation values *yT* = (*y*1, *y*2, ··· , *yn*) follows a probability (density) distribution function *f*(*yi*|*θi*), where *θi* = (*θ*1*i*, *θ*2*i*, *θ*3*i*, *θ*4*i*) = (*μi*, *σi*, *υi*, *τi*) is a parametric vector. The first two parameters of the model are usually defined as position and scale parameters, which are the mean vector and the mean variance vector of random variables. If there are other parameters in the distribution, such as the skewness vector *υi* and the kurtosis vector *τi* of random variable series, they are all designated shape parameters. In this paper, we only consider the first two parameters. A GAMLSS model can be expressed as a known monotonic link function which demonstrates the explanatory variables and random effects:

$$f(y\_i|\theta\_i) = X\_i \theta\_i \tag{1}$$

where *θi* are n-length vectors, and *θi* = (*θ*1*i*, *θ*2*i*, ··· , *θni*) *T*, *βi* = - *β*1*i*, *β*2*i*, ··· , *βIii T* is a parameter vector of length *Ii*, and *Xi* is an explanatory matrix of order *n* × *Ii*.

There are two basic algorithms for parameter estimation in the GAMLSS models. The first is the CG algorithm [25], and the other is the RS algorithm [9,10]. The latter algorithm is more suited for fitting mean and dispersion additive models, thus, it was used herein for parameter estimation.

In order to estimate the parameters in GAMLSS models, a penalized likelihood function *L* is usually introduced:

$$L = \sum\_{i=1}^{n} \log f(y\_i | \theta\_i) \tag{2}$$

Generally, the objective function is that the logarithmic likelihood function takes its maximum value. Then the regression parameter vector *βi* can be estimated by RS algorithm.

Analyzing the normality and independence of the residuals of the models can verify the quality of model fitting when there are no validated models. The model should be adequate with their mean nearly zero, variance nearly one, coefficient of skewness and kurtosis close to zero and three, respectively, and the Filliben coefficient [26] greater than the critical value given a certain sample size.

#### *3.2. Univariate Time Varying Model Based on GAMLSS Theory*

Four two-parameter distributions can reflect the actual hydrological process, and have thus been adopted within the framework of GAMLSS. The parametric model is given as:

$$f(y\_i|\mu) = X\_1\beta\_{1\prime} \; f(y\_i|\sigma) = X\_2\beta\_{2\prime} \tag{3}$$

*Water* **2018**, *10*, 819

The four two-parameter distributions for nonstationary flood frequency analysis are expressed as following:

(1) Gumbel distribution

$$f\_Y(y|\mu, \sigma) = \frac{1}{\sigma} \exp\left[\left(\frac{y-\mu}{\sigma}\right) - \exp\left(\frac{y-\mu}{\sigma}\right)\right], \quad -\infty < y < \infty, \ -\infty < \mu < \infty, \ \sigma > 0 \tag{4}$$

(2) Weibull distribution

$$f\_Y(y|\mu, \sigma) = \frac{\sigma y^{\sigma - 1}}{\mu^{\sigma}} \exp \left[ - \left( \frac{y}{\mu} \right)^{\sigma} \right], \ y > 0, \ \mu > 0, \ \sigma > 0 \tag{5}$$

(3) Gamma distribution

$$f\_Y(y|\mu, \sigma) = \frac{1}{(\sigma^2 \mu)^{1/\sigma^2}} \frac{y^{1/\sigma^2 - 1} e^{-y/(\sigma^2 \mu)}}{\Gamma(1/\sigma^2)}, \ y > 0, \ \mu > 0, \ \sigma > 0 \tag{6}$$

#### (4) Log-normal distribution

$$f\_Y(y|\mu, \sigma) = \frac{1}{\sqrt{2\pi\sigma^2}} \frac{1}{y} \exp\left\{-\frac{\left[\log(y) - \mu\right]^2}{2\sigma^2}\right\}, \ y > 0, \ \mu > 0, \ \sigma > 0. \tag{7}$$

#### *3.3. Bivariate-Joint Time Varying Model Based on Copulas*

For modeling the dependence structure of two or more random variables, copula functions are efficient mathematical tools. They were proposed first by Sklar [27], and have been widely used over the last decades. Considering a pair of random variables *X* and *Y*, with marginal distribution functions *u* = *FX*(*x*) = *P*(*X* ≤ *x*) and *v* = *FY*(*y*) = *P*(*Y* ≤ *y*), there will be a copula function *C* to describe the associated relationship, which can be expressed as:

$$F\_{X,Y}(x,y) = \mathbb{C}[F\_X(x), F\_Y(y)] = \mathbb{C}(u,v) \tag{8}$$

where *FX*,*Y*(*x*, *y*) is a joint cumulative distribution function (cdf) with margins *u* and *v*, all (*u*, *v*) ∈ (0, 1) 2 [28].

One kind of frequently-used copula is the Archimedean, which has three types, written as:

(1) Gumbel-Hougaard copula

$$\mathcal{C}\_{\theta}(u,v) = \exp\left\{-\left[\left(-\ln u\right)^{\theta} + \left(-\ln v\right)^{\theta}\right]^{\frac{1}{\theta}}\right\}, \; \theta \in [-1, \infty);\tag{9}$$

(2) Clayton copula

$$\mathbb{C}\_{\theta}(u,v) = \left(u^{-\theta} + v^{-\theta} - 1\right)^{\frac{1}{\theta}}, \theta \in (0, \infty);\tag{10}$$

(3) Frank copula

$$\mathbb{C}\_{\theta}(\boldsymbol{u}, \boldsymbol{v}) = -\frac{1}{\theta} \ln \left\{ 1 + \frac{[\exp(-\theta \boldsymbol{u}) - 1][\exp(-\theta \boldsymbol{v}) - 1]}{\exp(-\theta) - 1} \right\}, \theta \in \mathbb{R}. \tag{11}$$

The nonstationary models in this study were constructed through copulas, composed of two marginal distributions and a copula parameter *θ*. The marginal distributions were determined by the best nonstationary univariate models mentioned in Section 3.1, and only the copula parameter

*θ* needed to be solved herein. In literature the copula parameter *θ* was estimated by the inference function of margins method (IFM), which is given as:

$$L(u, v; \theta) = \sum \ln \mathbb{C}[F\_X(x), F\_Y(y); \theta]. \tag{12}$$

Let *∂L*/*∂θ* = 0, then *θ* can be calculated.

Two steps were carried out for model selection. Firstly, the Kolmogorov-Smirnov (K-S) test method was adopted to conduct the test of fit for copulae. Secondly, the models which passed K-S test were selected optimally according to the goodness-of-fit (GoF), which was represented by ordinary least squares (OLS) and Akaike information criterion (AIC) [29]. The K-S test statistic D and estimates of OLS and AIC were calculated by the following formulae:

$$D = \max\left\{ \left| \mathbf{C}\_k - \frac{i}{n} \right|, \left| \mathbf{C}\_k - \frac{i-1}{n} \right| \right\} \tag{13}$$

where *Ck* is the copula value of the measured sample series, *i* is the number of samples which meet the requirements that *x* ≤ *xi*, *y* ≤ *yi*, and *n* is the length of the sample series.

$$OLS = \sqrt{\frac{1}{n} \sum\_{i=1}^{n} (Pc\_i - P\_i)^2} \,\tag{14}$$

$$P\varepsilon\_i = \frac{i}{n+1'} \tag{15}$$

where *Pei* and *Pi* are the empirical frequency and theoretic frequency of measured sample series, respectively.

$$AIC = n \ln \left( \frac{1}{n} \sum\_{i=1}^{n} (Pe\_i P\_i)^2 \right) + 2m \tag{16}$$

where *m* is the number of model parameters. When the value of AIC is smaller, the model fitting is better.

The concept of return period in stationary frequency analysis is prone to misconceptions and misuses. New methods have been adopted to solve this problem, but have not worked well enough [30]. Since the return period is based on the probability, for the present study we explore the effect of nonstationarity on flood data focusing on the exceedance probability. As the return period is based on the probability, in the present work we studied the effect of nonstationarity on flood data based on exceedance probability. In the bivariate model, the OR-joint exceedance probability describes that at least one of the hydrologic variables *X* and *Y* exceed the values *x* and *y* respectively, that is, *P*∪ = *P*(*X* > *x* ∪ *Y* > *y*). The AND-joint exceedance probability describes that *X* and *Y* both exceed the values *x* and *y*, that is, *P*∩ = *P*(*X* > *x* ∩ *Y* > *y*). The *P*∪ and *P*∩ are written as:

$$P^{\downarrow} = 1 - \mathbb{C}[F\_X(\mathbf{x}), F\_Y(\mathbf{y})] \tag{17}$$

$$P^\cap = 1 - F\_X(\mathbf{x}) - F\_Y(\mathbf{y}) + \mathbb{C}[F\_X(\mathbf{x}), F\_Y(\mathbf{y})] \tag{18}$$

For a given data set, all the copula *C*(*u*, *v*) at the same probability level have the same exceedance probability. However, at least one combination of a given probability is more likely than others, namely the most-likely events. Therefore, the most-likely events can be selected as the point with the largest joint probability on the level curve, which was given by Gräler et al. [31]:

$$\begin{aligned} f(u,v) &= \underset{C\_{UV}(u,v) = k}{\text{argmax}} \, f\_{XY}(F\_X^{-1}(u), F\_Y^{-1}(v)) \\ &\ge = F\_X^{-1}(u), \, y = F\_Y^{-1}(v) \end{aligned} \tag{19}$$

where *k* is a given value of copula, and *x* and *y* can be calculated by the inverse cdf of the marginal distributions.

#### **4. Results**

In this study, following the identification of the nonstationarity for the Wangkuai Reservoir inflow flood sequences, two nonstationary models based on GAMLSS theory were established, in which the flood peak Q and flood volume W1 were considered as the independent response variables, and time t was adopted as the explanatory variable.

#### *4.1. Identification of Nonstationarity for Flood Sequences*

Firstly, we adopted the pettitt test and Mann-Kendall test to detect the change point and trend of the annual maximum flood peak Q and one-day flood volume W1 time series. Through pettitt test, the possible change points for the Q sequence are 1979 and 1996, and for W1 1979 and 1982. Since the test probability P of 1979 is the largest, and it is among the change points of both Q and W1 sequences, the most possible change point is 1979. This agrees with the results obtained by other researchers using different methods [32].

Before the trend analysis, the autocorrelation analysis of the flood sequences should be carried out. Since the autocorrelation coefficients of both the Q and W1 sequences are less than 0.1, it is considered that the autocorrelation of these sequences is not significant, so the trend test can be conducted directly. Without considering the change point, the non-parametric Mann-Kendall test was used to analyze the trend of the flood sequences. The statics *Un* of flood sequences Q and W1 are both less than −1.96, which shows that the flood sequences have passed the test at a significance of 5% and present a downward trend.

In order to consider the influence of the change point on the trend test, the Q and W1 sequences were divided into two subsequences and tested by non-parametric Mann-Kendall test method, respectively. Figure 2 shows the test results of the subsequences: although the subsequence of Q before 1979 shows a slight upward trend, and the subsequences of Q and W1 after 1979 show a downward trend, none of them passed the test of significance. Meanwhile, no trend has been found for the sequence of W1 before 1979.

In the case of full flood sequences presenting a significant downward trend, while the subsequences with 1979 as the dividing point showing no significant trend, the latter is more suitable for describing the inconsistency of the flood sequences.

**Figure 2.** Trend test results of the flood subsequences with consideration of change point.

#### *4.2. Nonstationary Model of Univariate Flood Frequency Analysis with Time as Covariate*

#### 4.2.1. Model Fitting Evaluation

A univariate nonstationary model with different distribution functions based on GAMLSS with time as covariate has been constructed based on the inflow flood sequences of Wangkuai Reservoir. In order to avoid an over-fitting problem due to excessive freedom degree, this study only considered the linear relationship between the distribution parameters and time *t*. The parameters include the distribution parameter *θ*1 (mean value of flood sequence) and *θ*2 (variance of flood sequence). The optimal fitting distribution, the optimal covariates of the distribution parameters, and the functional relationship between the distribution parameters and the optimal covariates were determined by the AIC criterion. Table 1 shows that, for flood peak (Q) and flood volume (W1) sequences, Weibull distribution and Gamma distribution had similar fitting effect, while log-normal distribution performed best with a minimum value of AIC. The functional relationships are written as Equation (20) for Q and Equation (21) for W1. For all the distribution models, the optimal covariates of the distribution parameters have been proved to satisfy the significant level 0.05 via *χ*2 test. Besides, for both Q and W1 time series, the distribution parameter *θ*1 presents a linear dependent relationship with time *t*, whereas *θ*2 is constant. Therefore, it can be concluded that the variation of the inflow flood sequences of Wangkuai Reservoir is mainly reflected by the mean value rather than the variance.

$$
\theta\_1 = 77.769 - 0.036t, \; \theta\_2 = 1.030 \text{ (for Q)}\tag{20}
$$

$$
\theta\_1 = 84.876 - 0.039t, \; \theta\_2 = 1.193 \, (for \, W\_1) \tag{21}
$$

Based on the above analysis, log-normal distribution was selected as the optimal distribution of both flood peak and volume. The corresponding residuals of the optimal distributions and the Filliben coefficients are shown in Table 2. For a sample size of 49, when the Filliben coefficients are both greater than 0.975, they satisfy the significant level 0.05 via *T* test. Thus, the model residuals of both Q and W1 sequences are acceptable and in normal distribution.

**Table 1.** The functional relationships of the explanatory variables and distribution parameters.

**Table 2.** The residuals of the optimal distributions and corresponding Filliben coefficients.

#### 4.2.2. Model Results Analysis

The summary of the associated results of the univariate nonstationary model is shown in Figure 3. Most of the observed data points are distributed in the grayscale range between 5% to 95% quantiles, indicating that the model results are able to capture the nonstationarity of the flood data. With time *t* as covariate, flood peak and volume time series show a downward trend over time. Especially when the quantile is larger (such as 95% quantile), the decrease trend is more significant. This proves

that variation of flood sequences has occurred under the changing environment, and the traditional "stationarity" hypothesis has become questionable. Thus, use of the traditional hydrological frequency analysis method may result in inaccurate results. The change of the flood sequences are mainly caused by climate change, land use change, and construction of Water conservancy projects [33–35]. Besides, the downward trend is significant before 1980, and tends to be gentle after 1980. This may be due to the change of land use. As revealed by Li et al. [34], water area, agricultural land, and grassland area decreased, while the forest area increased significantly in the control basin of Wangkuai Reservoir during 1970–1980, but from 1980 to 2000, the land use types kept almost invariant.

Figures 4 and 5 show the worm plot and normal QQ of the peak discharge Q and one-day volume W1. As can be seen, no significant departure from normality has been highlighted. Therefore, the results of residuals supported the inference that log-normal distribution provides a good fit to both Q and W1.

**Figure 3.** Summary of results of the univariate nonstationary model with Generalized Additive Models in Location, Scale, and Shape (GAMLSS) implementation. Red points represent the observed time series of peak discharge Q (**A**) and flood volume W1 (**B**). The solid black line represents 0.5 quantile; the dark grey region is the area between 0.25 and 0.75 quantiles; the light grey region is the area between 0.05 and 0.95.

**Figure 4.** Worm plot of the residuals of non-stationary time-varying model with GAMLSS implementation (Q in **A**, W1 in **B**). For a satisfactory fit of worm plot, all the observations should fall inside the two elliptic curves.

**Figure 5.** QQ plot of the residuals of non-stationary time-varying model with GAMLSS implementation (Q in **A**, W1 in **B**).

#### *4.3. Nonstationary Model of Bivariate-Joint Flood Frequency Analysis with Time as Covariate*

4.3.1. Parameters Calculation, Fitting Test, and Result Optimization of Joint Probability Distribution

With the above two log-normal distributions (Equations (19) and (20)) as marginal distributions, the joint distribution model of flood peak and volume based on copulas have been constructed according to Equations (9)–(11). The parameter of Copula functions *θ*, K-S test statistic D, and values of OLS and AIC were estimated by Equations (12)–(16). A significance level of *α* = 0.05 was used for the K-S test, and the corresponding standard quantile was approximately 0.1943 for *n* = 49. This means that when D-value is less than 0.1943, it passes the test. The results of the calculation, test, and optimization are presented in Table 3.

**Table 3.** The estimation of copula parameters and the associated fitting test and model selection.

It can be observed in Table 3 that Gumbel-Hougaard copula function and Clayton copula function have passed the K-S test, while the Frank Copula has not passed the test. Among them, the Gumbel-Hougaard copula function has the minimum values of OLS and AIC, indicating that the Gumbel-Hougaard copula function is more suitable to describe the extreme sequence of hydrological variables, such as the joint distribution of flood peak and volume. Therefore, the G-H Copula with parameter *θ* = 2.9224 was considered as the optimal copula for the bivariate-joint distribution of the Q and W1 time series of Wangkuai Reservoir.

4.3.2. Most Likely Combination of the Flood Peak and Volume

The "most likely" design events derived from Equation (11) have their maximum value of likelihood function on each probability-isoline, and all combinations of Q and W1 were illustrated in Figure 6. They have the similar undulating behavior as the univariate models shown in Figure 3.

**Figure 6.** Results of the most likely design events (combinations of flood peak and volume). Case illustration for *P*∪ = 0.02, 0.01 (**A**,**B**) and *P*∩ = 0.02, 0.01 (**C**,**D**).

Figure 6 shows that, under changing environmental conditions, with time *t* as the covariate of the marginal distribution, the combination values of the flood peak and volume under the same OR-joint and AND-joint exceedance probability both display a decreasing trend. Meanwhile, the combination values of the nonstationary model intersect with the traditional stationary results between 1970–1980. That is, before 1970, the most likely combination values considering the variation of distribution parameters over time were larger than fixed parameters (stationary), however, they were smaller after 1980.

#### **5. Discussion and Conclusions**

The changing environment makes the assumption of stationarity of flood sequences questionable. In this context, this paper constructed both univariate and bivariate nonstationary models with time as covariate based on GAMLSS theory for flood frequency analysis. The inflow flood peak Q and flood volume W1 series of Wangkuai Reservoir were used as basic data.

Within the framework of nonstationary flood frequency, this paper adopted four two-parameter distributions (Gumbel, Weibull, Gamma, and Log-Normal) as alternative distributions, which have the characteristics of power distribution and exponential distribution simultaneously. In the univariate nonstationary model with time as covariate, log-normal distribution performed best according to AIC criterion. The flood peak and volume time series presented a decreasing trend over time. Especially when the quantile is high (such as 95% quantile), the downward trend is more significant. Besides, the decreasing trend is significant before 1980, and tends to be gentle after 1980. This proves that variation of flood sequences has occurred under the changing environment.

Based on the optimal univariate models, copula functions were addressed to construct the dependence structure of Q and W1, with the two optimal log-normal distributions as marginal distributions. The results showed that only the Gumel-Hougaard copula can provide the best joint distribution. The most likely events have similar undulating behavior to the univariate models, and the combination values of the flood peak and volume under the same OR-joint and AND-joint exceedance probability both display a decreasing trend. Meanwhile, the combination values of the nonstationary model intersect with the traditional stationary results between 1970–1980. That is, before 1970, the most likely combinations considering the variation of distribution parameters over time were larger than fixed parameters (stationary), whereas it became the opposite after 1980.

**Author Contributions:** Conceptualization, P.F.; Methodology, T.Z.; Software, S.T.; Validation, Y.W.; Formal Analysis, T.Z.; Investigation, B.W.; Resources, B.W.; Data Curation, Y.W.; Writing Original Draft Preparation, T.Z.; Writing Review & Editing, T.Z.; Visualization, S.T.; Supervision, P.F.; Project Administration, P.F.; Funding Acquisition, T.Z.

**Funding:** This research was funded by National Natural Science Foundation of China grant number 51609165, China Postdoctoral Science Foundation grant number 2017T100159, and the Foundation of State Key Laboratory of Hydraulic Engineering Simulation and Safety (Tianjin University) grant number HESS-1607.

**Acknowledgments:** We are grateful to Hydrology and Water Resource Survey Bureau of Hebei Province for providing the hydrometeorological data. The authors also acknowledge the reviewers for providing valuable advices to improve our papers.

**Conflicts of Interest:** The authors declare no conflict of interest. The funding sponsors had no role in the design of the study; in the collection, analyses, or interpretation of data; in the writing of the manuscript, and in the decision to publish the results.

#### **References**

© 2018 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (http://creativecommons.org/licenses/by/4.0/).


**4. Search-based replanning on time-varying environment**

In real-world application, there is often a scenario that the robot initially does not know a priori information about its surroundings. We cannot encode the world space information each time the robot runs, because it is expensive, tedious, and infeasible due to rapid changes in practice. To maintain collision-free path, one can naively rerun A\* to replan the shortest path from the point that the robot detects changes. However, this naïve approach will waste computation by reprocessing cells that are irrelevant to compute a new path and hence increase idle time between each search. This section will demonstrate search-based algorithms to solve mentioned problem in time-variant environment.


**Conclusion**

IQoro training significantly improves all the symptoms of hiatus hernia, likely through improved hiatal competence.

*<IP_ADDRESS> Study: Oral neuromuscular training relieves hernia-related dysphagia and GERD symptoms as effectively in obese as in non-obese patients*


**2. Materials and Methods**

The identification and determination of heat waves over the study area (Greece) was achieved using the computation of the excess heat factor (EHF) index, which provides a measure of the environmental temperature load [10] and the intensity of a potential heat wave [12]. The EHF is the product of the multiplication of two other excess heat indices (EHIs), namely, EHIaccl and EHIsig.

$$\text{EHF} = \text{EHI}\_{\text{sig}} \times \max\{1, \text{EHI}\_{\text{acccl}}\} \tag{1}$$

EHIaccl is defined as:

$$EHIaclcl = \left[\frac{T\_i + T\_{i-1} + T\_{i-2}}{3}\right] - \left[\frac{T\_{i-3} + \dots \, T\_{i-32}}{30}\right],\tag{2}$$

where the first term of Equation (2) is the average daily Tmean for a three-day period and the second term is the average daily Tmean of the preceding thirty days. As proposed by Nairn et al. [34], this is an acclimatization index and its positive (negative) values are related to hot (cold) weather conditions. It determines a period of heat that is warmer than the recent past [10], and it should be highlighted that this index is not influenced by the potential general warming trend [34]. This index describes an important factor of the influence of heat to the population because, even though humans tend to acclimatize themselves to their environmental local climate according to the temperature variations throughout the year, they may be unprepared to an abrupt temperature rise above that of the recent past [12]. Thus, positive EHIaccl values indicate a lack of acclimatization to the warmer temperatures which may result in negative health impacts.

The second term of the EHF equation is the significance index (EHIsig) defined by the following equation:

$$EHIsig = \left[\frac{T\_i + T\_{i-1} + T\_{i-2}}{3}\right] - T\_{\mathfrak{F}\mathfrak{F}\_{\prime}}\tag{3}$$

which is calculated by the difference of the average daily Tmean for a three-day period minus the 95th percentile of the daily Tmean. The percentile is computed for a reference period of 30 years (1971–2000) using the daily values of the mean temperatures for all days throughout the year. A heat wave occurs when EHIsig is positive, while the comparison with the 95th percentile measures the statistical significance of the heat event [30]. The authors also underline the fact that, since T95 is computed for a fixed climatological period, the EHIsig (contrary to the EHIaccl) is expected to become more extreme under a general warming trend.

Therefore, the excess heat factor expresses the long-term temperature anomalies amplified by the short-term ones [10], and the days with positive EHF values indicate heat wave conditions, while the higher the values of the index, the more intense is the heat wave. However, according to Perking and Alexander [1], a heat wave episode will be defined when at least three consecutive days present EHF values above zero.

Finally, daily Tmean in all the above indices should be computed by averaging the Tmin and Tmax daily values since the diurnal temperature variation is highly associated with the ability of the biological systems to recover from high heat load. Hence, for this study, the daily Tmin and Tmax time series, derived from 14 meteorological stations, were used for a 55-year period starting from 1958 to 2012, for the computation of the Tmean values. Except for the data for the Thessaloniki station which were available from the meteorological station of the Aristotle University of Thessaloniki (AUTh), the remaining station data were provided by the Hellenic National Meteorological Service (HNMS). These data were proven to be homogenous according to the Alexandersson test [38] and had no gaps. The geographical distribution of the station locations is presented in Figure 1.

It is worth mentioning at this point that although abnormally warm conditions may occur even during the winter months, it was decided in this case to compute the EHF index for the hottest period of the year, that is, June, July, August, and September (the JJAS period) since summer heat waves tend to be more intense with severe impacts for humans during these months. Finally, adopting the definition used by Perkins and Alexander [1] who mentioned that a heat wave occurs when the abnormally hot conditions (EHF > 0) persist for more than three consecutive days, the duration of spells longer than 3 days was calculated for the 14 stations under study during the 55-year time period.

**Figure 1.** Geographical distribution of the stations under study. For each station, the altitude where it is located can be found on the right of the map. For each station, the average Tmin, Tmax, and Tmean values for the period 1958–2012 are provided. The sign of the trend of these time series is found in the brackets (). The asterisk indicates the statistical significance of the trends at the 95% level of significance.


**1. Introduction**

170 Hepatocellular Carcinoma – Basic Research

[18] Konduri SD, Srivenugopal KS, Yanamandra N*, et al*: Promoter methylation and silencing of

[20] Hube F, Reverdiau P, Iochmann S, Rollin J, Cherpi-Antar C and Gruel Y:

[21] Rollin J, Iochmann S, Blechet C*, et al*: Expression and methylation status of tissue factor pathway inhibitor-2 gene in non-small-cell lung cancer. Br J Cancer 92: 775-783, 2005. [22] Guo H, Lin Y, Zhang H*, et al*: Tissue factor pathway inhibitor-2 was repressed by CpG

[23] Nobeyama Y, Okochi-Takada E, Furuta J*, et al*: Silencing of tissue factor pathway inhibitor-2 gene in malignant melanomas. Int J Cancer 121: 301-307, 2007. [24] Wong CM, Ng YL, Lee JM*, et al*: Tissue factor pathway inhibitor-2 as a frequently silenced tumor suppressor gene in hepatocellular carcinoma. Hepatology 45: 1129-1138, 2007. [25] Kempaiah P, Chand HS and Kisiel W: Identification of a human TFPI-2 splice variant that is upregulated in human tumor tissues. Mol Cancer 6: 20, 2007. [26] Tasiou A, Konduri SD, Yanamandra N*, et al*: A novel role of tissue factor pathway

metalloproteinases in human glioma cells. Oncogene 22: 4509-4516, 2003. [19] Pulukuri SM, Gorantla B and Rao JS: Inhibition of histone deacetylase activity promotes

activator. J Biol Chem 282: 35594-35603, 2007.

cancer cells. BMC Mol Biol 8: 110, 2007.

Metastasis 18: 303-308, 2000.

vitro. Int J Oncol 18: 127-131, 2001.

and in vivo. Gynecol Oncol 83: 325-333, 2001.

meningioma cell line. Int J Oncol 29: 25-32, 2006.

kDa MSPI. Int J Cancer 76: 749-756, 1998.

Invest 107: 1117-1126, 2001.

choriocarcinoma cells. Biol Chem 384: 1029-1034, 2003.

the tissue factor pathway inhibitor-2 (TFPI-2), a gene encoding an inhibitor of matrix

invasion of human cancer cells through activation of urokinase plasminogen

Transcriptional silencing of the TFPI-2 gene by promoter hypermethylation in

hypermethylation through inhibition of KLF6 binding in highly invasive breast

inhibitor-2 in apoptosis of malignant human gliomas. Int J Oncol 19: 591-597, 2001. [27] Konduri SD, Tasiou A, Chandrasekar N, Nicolson GL and Rao JS: Role of tissue factor

pathway inhibitor-2 (TFPI-2) in amelanotic melanoma (C-32) invasion. Clin Exp

pathway inhibitor-2 (TFPI-2), decreases the invasiveness of prostate cancer cells in

2/MSPI decreases the invasive potential of human choriocarcinoma cells in vitro

(rAAV) expressing TFPI-2 inhibits invasion, angiogenesis and tumor growth in a

inhibits invasion and tumor growth in vitro and in vivo in a malignant

invasion are inhibited by the matrix-associated serine protease inhibitor TFPI-2/33

activation by tissue factor pathway inhibitor-2/matrix-associated serine protease

the production of active matrix metalloproteinase-2 and is down-regulated in cells

inhibitor of matrix metalloproteinases with implications for atherosclerosis. J Clin

[28] Konduri SD, Tasiou A, Chandrasekar N and Rao JS: Overexpression of tissue factor

[29] Jin M, Udagawa K, Miyagi E*, et al*: Expression of serine proteinase inhibitor PP5/TFPI-

[30] Yanamandra N, Kondraganti S, Gondi CS*, et al*: Recombinant adeno-associated virus

[31] Kondraganti S, Gondi CS, Gujrati M*, et al*: Restoration of tissue factor pathway inhibitor

[32] Rao CN, Cook B, Liu Y*, et al*: HT-1080 fibrosarcoma cell matrix degradation and

[33] Rao CN, Mohanam S, Puppala A and Rao JS: Regulation of ProMMP-1 and ProMMP-3

[34] Izumi H, Takahashi C, Oh J and Noda M: Tissue factor pathway inhibitor-2 suppresses

[35] Herman MP, Sukhova GK, Kisiel W*, et al*: Tissue factor pathway inhibitor-2 is a novel

human glioblastoma cell line. Int J Cancer 115: 998-1005, 2005.

inhibitor. Biochem Biophys Res Commun 255: 94-98, 1999.

harboring activated ras oncogenes. FEBS Lett 481: 31-36, 2000.

The last four decades witnessed a brilliant career of proteoglycans (PGs). Once regarded as mere space-fillers or passive structural components of matrices and charge-selective barriers, these fascinating molecules have been increasingly acknowledged as key players in cell-cell and cell-matrix communication, and have become recognized as modulators of most, if not all, aspects of cell behavior including survival, proliferation, and migration. Simultaneously, the range of disease processes with known involvement of PGs has steadily expanded, now covering areas as diverse as host-pathogen interactions, regulation of pathologic fibrogenesis, and tumor progression. Characteristic alterations of PGs in various human malignant tumors, including HCC, were first described more than 20 years ago (for an early review, see Tímár & Kovalszky, 1995).

PGs, glycanated proteins with extensive posttranslational modifications, consist of a protein core and one or more long, linear, sulfated polysaccharide chains, called glycosaminoglycans (GAGs). GAGs are ligated to the protein core at specific serine, threonine, or asparagine residues, although the exact signal sequences that designate the position of attachment are mostly unknown. The multifunctionality of PGs arises from their inherently complex structure: some functions are assigned to the core protein, while others are fulfilled by the GAG chains.

The synthesis of each GAG chain (recently reviewed by Ly et al., 2010) is introduced by the attachment of a short linkage region to Ser in the case of heparan sulfate (HS)/heparin and chondroitin sulfate/dermatan sulfate (CS/DS), and either Asn or Ser/Thr in the case of keratan sulfate (KS) type I and type II, respectively. During the elongation phase of GAG synthesis, acetylated hexosamine and hexuronic acid or galactose residues are added in an alternating fashion to the growing polysaccharide chain. GAGs are classified by their disaccharide composition: the dimeric building block is N-acetyl-glucosamine / glucoronic acid in HS and heparin; N-acetyl-galactosamine / glucoronic acid in CS and DS; and Nacetlyl-glucosamine / galactose in KS. Completed GAG chains then undergo various chemical modifications including *N*-deacetylation, *N*- and *O*-sulfation, and epimerization of the hexuronic acid. Heparin, for example, differs from HS in the extent of sulfation (heparin is sulfated uniformly and nearly exhaustively, whereas HS is sulfated only partially and in a

Proteoglycans in Chronic Liver Disease and Hepatocellular Carcinoma: An Update 173

patterned manner); and DS differs from CS in the degree of uronic acid epimerization (0% in CS vs. 1-100% in DS.) The extent and pattern of modifications not only vary between different GAGs and PGs, but also depend on the type and actual state of a cell, which contributes a great deal to the biological diversity of PGs. If this were still not enough of versatility, PGs may undergo further editing once they are in place: in the matrix or on the cell surface, they may be subject to the action of endoglycosidases that cleave the GAG chain, proteases that cut the protein core, and endosulfatases capable of removing sulfate groups from internal sugar residues.

Historically, PGs were sorted by the type of their GAG chain into one of the categories HSPG, CS/DSPG, or KSPG. Later, however, the discovery that several PGs carry more than one type of GAG (i.e., syndecans and betaglycan carry both HS and CS; aggrecan carries both CS and KS II) prompted a new classification based on structure and tissue localization. In this revised system, each PG belongs to one of three major families: 1) small leucine-rich proteoglycans or SLRPs; 2) modular PGs, further divided into a) hyalectans or hyaluronan-binding PGs and b) non-hyaluronan-binding PGs of the basement membrane; and 3) cell surface PGs. Nevertheless, both the old and new classifications fall short of being perfect; neither is free of overlaps, and neither can properly accommodate, for example, serglycin or endocan. In this review, we shall follow a sort of "hybrid" classification that fits best for our purposes.

A complete listing of all currently known PGs seems unnecessary here (for a comprehensive review, the Reader is referred to Esko et al., 2009); this paper is restricted in scope to PGs present in the healthy or diseased liver, and will concentrate on those involved in, or affected by, chronic liver disease and hepatocarcinogenesis. Also, with a focus on human disease, PGs reported to be present in the liver of experimental animals but not of humans will be omitted.


*Article* **Dynamic Transcriptome Analysis of Anther Response to Heat Stress during Anthesis in Thermotolerant Rice (***Oryza sativa* **L.)**

**Gang Liu** †**, Zhongping Zha** †**, Haiya Cai, Dandan Qin, Haitao Jia, Changyan Liu, Dongfeng Qiu, Zaijun Zhang, Zhenghuang Wan, Yuanyuan Yang, Bingliang Wan, Aiqing You and Chunhai Jiao \***

Hubei Key Laboratory of Food Crop Germplasm and Genetic Improvement, Food Crops Institute, Hubei Academy of Agricultural Sciences, Wuhan 430064, China<EMAIL_ADDRESS>(G.L<EMAIL_ADDRESS>(Z.Z<EMAIL_ADDRESS>(H.C<EMAIL_ADDRESS>(D.Q<EMAIL_ADDRESS>(H.J<EMAIL_ADDRESS>(C.L<EMAIL_ADDRESS>(D.Q<EMAIL_ADDRESS>(Z.Z<EMAIL_ADDRESS>(Z.W.);

<EMAIL_ADDRESS>(Y.Y<EMAIL_ADDRESS>(B.W<EMAIL_ADDRESS>(A.Y.)

Received: 13 January 2020; Accepted: 7 February 2020; Published: 10 February 2020

**Abstract:** High temperature at anthesis is one of the most serious stress factors for rice (*Oryza sativa* L.) production, causing irreversible yield losses and reduces grain quality. Illustration of thermotolerance mechanism is of great importance to accelerate rice breeding aimed at thermotolerance improvement. Here, we identified a new thermotolerant germplasm, SDWG005. Microscopical analysis found that stable anther structure of SDWG005 under stress may contribute to its thermotolerance. Dynamic transcriptomic analysis totally identified 3559 differentially expressed genes (DEGs) in SDWG005 anthers at anthesis under heat treatments, including 477, 869, 2335, and 2210 for 1, 2, 6, and 12 h, respectively; however, only 131 were regulated across all four-time-points. The DEGs were divided into nine clusters according to their expressions in these heat treatments. Further analysis indicated that some main gene categories involved in heat-response of SDWG005 anthers, such as transcription factors, nucleic acid and protein metabolisms related genes, etc. Comparison with previous studies indicates that a core gene-set may exist for thermotolerance mechanism. Expression and polymorphic analysis of agmatine-coumarin-acyltransferase gene *OsACT* in different accessions suggested that it may involve in SDWG005 thermotolerance. This study improves our understanding of thermotolerance mechanisms in rice anthers during anthesis, and also lays foundation for breeding thermotolerant varieties via molecular breeding.

**Keywords:** rice; heat stress; transcriptome; anther; anthesis


**8. The Maslow's hierarchy angle**

The Maslow's hierarchy of needs is a motivational theory in psychology comprising a five- tier model of human needs, often depicted as hierarchical levels within a pyramid. Needs further south in the hierarchy must be satisfied before individuals can sort out needs northwards. From bottom-up, the needs are: physiological, safety, love and belonging, esteem and self-actualization.

Maslow [9] using the above hierarchy of human need, stated that people are motivated to achieve certain needs and that some needs take precedence over others. From the list, he postulated that our most basic need is that of physical survival and this will be the first factor that motivates our behavior. Articulating that money is key to solving our physical needs, it therefore means that man gravitates quickly to wherever he can acquire funds with ease. Understanding the preceding fact, whatever positions have great monetary benefits attached usually will attract vehement but detrimental interest. Political Leadership is usually not an exception. The monetary emoluments attached to political positions; particularly in developing countries is enormous. As a result, the race to political offices are highly contested to critical levels such that it usually becomes a very risky affair. The acute poverty level in such countries ironically, is connected to the political leaders abysmal performance, creates an intense pursuit for political positions as an easy route to economic hardship.

#### **9. Comparative analysis of emoluments in developing countries**

To state that the pay packages of political office holders completely dwarfs that of civil servants/workers in developing countries is not an overstatement. From

critical analysis and information sourcing, the civil servants remuneration is match for the average political leader in most countries in Africa. Citing scenarios from Nigeria, according to information captured from online tabloid-Pulse, the average Nigerian Senator earns a minimum of N32million monthly as pay package.  If we compare this with the earnings of the highest level researcher from the knowledge setup- a Professor, you will see the clear imbalance. The highest earning Professor in a country like Nigeria, with all emoluments put together earns nothing above the realms of N8million annually. Considering that they are the hope for the dissemination of new knowledge and discoveries as well as inventions, such pay in comparison with their political counterparts is appalling. This type of pay disparity will not encourage the youths and upcoming professionals to develop interest in other areas as much as politics. This what will bring all sorts on individuals into the political terrain. Rather than come for service, they will come in for wealth gain which never makes for productive leadership.


**1. Introduction**

50 Chromatography and Its Applications

Subbiah, V. (1999). *Method Of Isolating Cucurbitacin*.United States Patent. 1999, patent

Sur, S.V. (1991). Isolation of Oxygen-Containing Monoterpenoids of Essential Oils by

Svedström, U.; Vuorela, H.; Kostiainen, R.; Tuominen, J.; Kokkonen, J.; Rauha, J.P.; Laakso, I.

Syu, W.J.; Shen, C.C.; Don, M.J.; Ou, J.C.; Lee, G.H. & Sun, C.M. (1998). Cytotoxicity of

Tiuman, T.S.; Ueda-Nakamura, T.; Garcia Cortez, D.A.; Dias Filho, B.D.; Morgado-Diaz, J.A.;

Umlauf, D.; Zapp, J.; Becker, H. & Adam, K.P. (2004). Biosynthesis of the Irregular

Vortherms, T.A. & Roth BL. Salvinorin A: From natural Product to Human Therapeutics.

Wall, M.E. & Wani, M.C. (1996). Camptothecin and Taxol: From Discovery to Clinic. *Journal*

Wang, G., Tang, W. & Bidigare, R.R. In Ed: Zhang L., & Demain, A.L., (2005). *Natural*

Yamada, K.; Murata, T.; Kobayashi, K.; Miyase, T.; & Yoshizaki, F. (2010). A Lipase Inhibitor

Crataegus Leaves and Flowers. *Phytochemistry*, Vol. 60, pp. 821-825.

Preparative Adsorbtion Chromatography with Gradient Elution. *Khimiya*

& Hiltunen, R. (2002). Isolation and Identification of Oligomeric Procyanidins from

Curcuminoids and Some Novel Compounds from Curcuma zedoaria. *Journal of*

Souza W. & Nakamura, C.V. (2005). Antileishmanial Activity of Parthenolide, a Sesquiterpene Lactone Isolated from Tanacetum parthenium*. Antimicrobial Agents*

Monoterpene Artemisia Ketone, the Sesquiterpene Germacrene D and Other Isoprenoids in Tanacetum vulgare L. (Asteraceae). *Phytochemistry*, Vol. 65, pp. 2463-

*Products Drug Discovery and Therapeutic Medicine Terpenoids As Therapeutic Drugs As Pharmaceutical Agents.* Humana Press, New Jersey, United States of Amerika. Watanabe, K.; Mimaki, Y.; Sakagami, H.; & Sashida, Y. (2002). Cycloartane Glycosides from

the Rhizomes of Cimicifuga racemosa and Their Cytotoxic Activities. *Chemical*

Monoterpene and Monoterpene Glycosides from Monarda punctata.

number: 5,925,356.

2470.

*Prirodnykh Soedinenii*, Vol. 5, pp. 634-637.

*Natural Products*, Vol. 61, pp. 1531-1534.

*and Chemotherapy*, Vol. 49, pp. 176-182.

*Molecular Interventions* 2006 Oct;6(5):257-65.

*of Ethnopharmacology*, Vol. 51, pp. 239-254.

*Pharmaceutical Bulletin*, Vol. 50, pp. 121-125.

*Phytochemistry*, Vol. 71, pp. 1884-1891.

Herbal medicines are major source of aphrodisiacs and have been used worldwide for thousands of years by different cultures and civilizations. Recently, consumption of dietary supplements has been becoming more popular around the world. Unfortunately, the adulteration of dietary supplements with undeclared synthetic chemical compounds is steadily increasing according to the literature. Some herbal products advertised as "all natural" have in contrast been found to contain synthetic PDE-5 inhibitors. There are currently three PDE5 inhibitors **Sildenafil** (Langtry & Markham, 1999) (Viagra; Pfizer, New York, US), **Tadalafil** (Meuleman, 2003) (Cialis; Eli Lilly, Indianapolis, US), and **Vardenafil** (Keating & Scott, 2003) (Levitra; Bayer Pharmaceuticals Co, Wuppertal, Germany), approved worldwide for the treatment of male erectile dysfunction, further two agents **Udenafil** (Salem et al., 2006) (Zydena; Dong-A PharmTech Co, Korean), **Mirodenafil** (Jung, 2008) (Mvix, Life Science R&D Center of SK chemical, Beijing, Tianjin, Shanghai) were licensed only in Korea. They produce vascular smooth muscle relaxation, promote penile blood flow, and hence, induce erection. These kinds of commercially available herbal aphrodisiac products have been spiked with the above-mentioned legal drugs, but also with their analogues, which have not been subjected to formal pharmacokinetic or other pharmacological testing in either humans or animals.

The practice of self-medication by an increasing number of patients, the incessant aggressive advertising of these herbal aphrodisiacs, the invasion of the medicinal market with uncontrolled dietary supplements and the absence of real directives amplifies the potential health hazards to the community. Since the sildenafil is an chemical, it must not been found in any foodstuffs, but an increasing number of sildenafil analogues have been discovered in dietary and herbal supplements even in soft drinks, this number is steadily increasing and some time their types are unknown and necrosis not observed instead of paracetamol. Hence, it is prudent to test the safety and efficacy thus might have unknown and harmful side-effects. Structural analogues are also synthetic chemicals with slightly altered chemical structures and have similar erectile effects on the body. Nevertheless, it is not uncommon

\* Corresponding Author

Chromatographic Separation and Identification of Sildenafil

Cl

HN

N N

Fig. 3. (a) Acetylvardenafil; (b) Desulfovardenafil; (c) Gendenafil.

O

O

minutes or less) (Bell & Palmer, 2011) Fig 2.

Fig. 2. Avanafil.

N

et al., 2007; Lai et al., 2007b) Fig 4(a).

2011) in Fig 4(c).

N

and Yohimbine Analogues Illegally Added in Herbal Supplements 53

**Avanafil** was discovered by pharmaceutical developer "Vivus" as PDE5 inhibitors, the safety and efficacy for erectile dysfunction, a possible lunch in 2012. It has a particular advantage over its potential competitors, the effects are demonstrable very quickly (in 15

N N

H

N N

O

(a) (b) (c)

Another new analogue of vardenafil, **Desulfovardenafil**, in which the *N*-ethylpiperazine ring and the sulphonyl group were removed from the vardenafil structure, was identified for enhancing erectile function in herbal health product marketed, namely Power58 Platinum (Lai et al., 2007, Lam, 2007) Fig 3(b). New sildenafil analogues **Gendenafil** had an acetyl group instead of sulfonyl-*N*-methylpiperazine moiety was determined as 5-[2-ethoxy-5-acetyl-phenyl]-1-methyl-3-*n*-propyl-1,6-dihydro-*7H*-pyrazolo[4,3-*d*]pyrimidin-7-one (Lin et al., 2008) Fig3(c). Other vardenafil analogue was found to be added illegally into a dietary supplement marketed for male erectile dysfunction (MED). Its structure was determined as 2-(2-ethoxyphenyl)-5-methyl-7-propyl-imidazo[5,1-*f*][1,2,4]triazin-4(*3H*)-one and called "**Piperidenafil** or **Piperidino vardenafil** or trivial name **Pseudovardenafil**" (Park

Very recently, two new analogues of sildenafil in which the piperazine ring and the sulfonyl group were replaced by a piperazinone and a hydroxyethyl structure, respectively were isolated from a herbal product in Germany. Based on the piperazinone structure, the compounds were named **Piperazinonafil** in Fig 4(b) and **Isopiperazinonafil** (Wollein et al.,

O

N H

The one of the last discovered illegal PDE-5 inhibitor is called as "**Acetylvardenafil**" which was found in dieatery supplement known as MEGATON**®** in USA (Lee et al., 2011).

O O N

Sulfonyl group of Vardenafil was substituted by an acetyl group Fig 3(a) .

N

O HN

OH

N

N

N

O

HN

N

O

O

N N

N

for chemicals with similar structures to possess slightly or entirely different before any new chemical is licensed as drug for human use. This testing process is lengthy and costly; on average, it takes 9.5 years and costs US\$800 million to license a new drug. Many drug analogues, without the aforementioned drug testing process, are available for human consumption properties. Phenacetin, structurally similar to paracetamol, has been associated with renal papillary necrosis not observed with paracetamol (Poon, 2007). Many drug analogues, without the aforementioned drug testing process, are avaliable for human comsuption via different channels. Examples include analogues of psychoactive drugs, anabolic steroids, and Sibutramine which was one of the most abused compound as antiobesity drugs.

The most commonly reported side effects of sildenafil are headaches, flushing of the face, upset stomach and nasal congestion. We met a dieatery supplement having combination of sildenafil and paracetamol, probably in order to prevent headeache caused by sildenafil. Other side effects include sensitivity to light, blurred vision, urinary tract infection, diarrhea and dizziness. The main problem with sildenafil and analogues are that they interact with many other medications. They can rapidly decrease blood pressure up to dangerously low. They can interact with nitrates such as nitrogliserin, which are often prescribed to heart patients. Sildenafil has not to be administrated to patients with heart problems taking nitrate medications because of the severe potentiation of vasodilatory effects.

It is well known, the first developed and consequently the most famous phosphodiesterase inhibitor is Sildenafil (Langtry & Markham, 1999) approved by the FDA in early April 1998. Novel PDE5 inhibitor, **Lodenafil carbonate,** breaks down in the body to form two molecules of the active drug lodenafil. This formulation has higher oral bioavailability than the parent drug (Toque et al., 2008) Fig 1.

Lodenafil

Fig. 1. Prodrug lodenafil carbonate.

**Avanafil** was discovered by pharmaceutical developer "Vivus" as PDE5 inhibitors, the safety and efficacy for erectile dysfunction, a possible lunch in 2012. It has a particular advantage over its potential competitors, the effects are demonstrable very quickly (in 15 minutes or less) (Bell & Palmer, 2011) Fig 2.

Fig. 2. Avanafil.

52 Chromatography and Its Applications

for chemicals with similar structures to possess slightly or entirely different before any new chemical is licensed as drug for human use. This testing process is lengthy and costly; on average, it takes 9.5 years and costs US\$800 million to license a new drug. Many drug analogues, without the aforementioned drug testing process, are available for human consumption properties. Phenacetin, structurally similar to paracetamol, has been associated with renal papillary necrosis not observed with paracetamol (Poon, 2007). Many drug analogues, without the aforementioned drug testing process, are avaliable for human comsuption via different channels. Examples include analogues of psychoactive drugs, anabolic steroids, and Sibutramine which was one of the most abused compound as anti-

The most commonly reported side effects of sildenafil are headaches, flushing of the face, upset stomach and nasal congestion. We met a dieatery supplement having combination of sildenafil and paracetamol, probably in order to prevent headeache caused by sildenafil. Other side effects include sensitivity to light, blurred vision, urinary tract infection, diarrhea and dizziness. The main problem with sildenafil and analogues are that they interact with many other medications. They can rapidly decrease blood pressure up to dangerously low. They can interact with nitrates such as nitrogliserin, which are often prescribed to heart patients. Sildenafil has not to be administrated to patients with heart problems taking nitrate

It is well known, the first developed and consequently the most famous phosphodiesterase inhibitor is Sildenafil (Langtry & Markham, 1999) approved by the FDA in early April 1998. Novel PDE5 inhibitor, **Lodenafil carbonate,** breaks down in the body to form two molecules of the active drug lodenafil.

 Prodrug lodenafil carbonate.

N S O O The one of the last discovered illegal PDE-5 inhibitor is called as "**Acetylvardenafil**" which was found in dieatery supplement known as MEGATON**®** in USA (Lee et al., 2011). Sulfonyl group of Vardenafil was substituted by an acetyl group Fig 3(a) .

Fig. 3. (a) Acetylvardenafil; (b) Desulfovardenafil; (c) Gendenafil.

Another new analogue of vardenafil, **Desulfovardenafil**, in which the *N*-ethylpiperazine ring and the sulphonyl group were removed from the vardenafil structure, was identified for enhancing erectile function in herbal health product marketed, namely Power58 Platinum (Lai et al., 2007, Lam, 2007) Fig 3(b). New sildenafil analogues **Gendenafil** had an acetyl group instead of sulfonyl-*N*-methylpiperazine moiety was determined as 5-[2-ethoxy-5-acetyl-phenyl]-1-methyl-3-*n*-propyl-1,6-dihydro-*7H*-pyrazolo[4,3-*d*]pyrimidin-7-one (Lin et al., 2008) Fig3(c). Other vardenafil analogue was found to be added illegally into a dietary supplement marketed for male erectile dysfunction (MED). Its structure was determined as 2-(2-ethoxyphenyl)-5-methyl-7-propyl-imidazo[5,1-*f*][1,2,4]triazin-4(*3H*)-one and called "**Piperidenafil** or **Piperidino vardenafil** or trivial name **Pseudovardenafil**" (Park et al., 2007; Lai et al., 2007b) Fig 4(a).

Very recently, two new analogues of sildenafil in which the piperazine ring and the sulfonyl group were replaced by a piperazinone and a hydroxyethyl structure, respectively were isolated from a herbal product in Germany. Based on the piperazinone structure, the compounds were named **Piperazinonafil** in Fig 4(b) and **Isopiperazinonafil** (Wollein et al., 2011) in Fig 4(c).

Chromatographic Separation and Identification of Sildenafil

N

Fig. 6. (a) *N*-octylnor-tadalafil; (b) Cyclopentynafil.

N H N

O

Fig. 7. Condensation product of aminotadalafil and hydroxymethylfuraldehyde.

N

N O

OH

O

O

O

Last flash development for the treatment of MED is discovering of **Zoraxel** (RX-10100) by Rexahn Pharmaceutical company (Albersen et al., 2010). Zoraxel is containing clavulanic acid that is centrally acting in the CNS and may be a more effective MED treatment for patients who are responsive or unresponsive to PDE-5 inhibitors. It is being developed as an orally administered, on-demand tablet to treat sexual dysfunction, and has extensive and well-established safety in humans. For the future, it is being expected, Zoraxel will be on to

O

O

O

N H

al., 2010) Fig 7.

worldwide best-selling drug.

N

O

and Yohimbine Analogues Illegally Added in Herbal Supplements 55

In Taiwan, one of the dietary supplement which was claimed on the treatment of male erectile dysfunction was firstly screened in 2009 and Tadalafil and its doctored version was newly identified. Since it is having amino group instead of methyl in tadalafil it was named as **aminotadalafil** (Zou et al., 2006; Lin et al., 2009) Table-1. This compound has two asymmetric carbons, theoritically two pairs of enantiomers exist. The chromatographic separation of its stereoisomers was reported by using chiral LC-MS (Kurita et al., 2008). Using this method, *RR*-Aminotadalafil and *SR*-Aminotadalafil were detected in some health food. In addition, an interaction product of aminotadalafil was isolated from an illegal health food product. The structure of the interaction product was elucidated and unknown compound was characterized as condensation product of aminotadalafil and hydroxymethylfuraldehyde and is probably the result of a drug-excipient incompatibility (Häberli et

HN

S N

N

(a) (b)

O O

N

O

O

N N

Fig. 4. (a) Piperidenafil; (b) Piperazinonafil; (c) Isopiperazinonafil.

Another sildenafil analogue was detected from a health supplement claimed for human MED, the structure of this new analogue was characterized as **dithiodesmethylcarbodenafil** containing 2 thiocarbonyl groups instead of 2 carbonyl groups, and 4-methyl substitution on the piperazine ring, rather than 4-ethyl substitution when compared to sildenafil (Ge et al., 2011) Fig. 5(a).

Fig. 5. (a) Dithio-desmethylcarbodenafil; (b) **N**itroso-prodenafil.

Other new unapproved analogue of sildenafil was detected in capsules of a herbal dietary supplement promoted as a libido enhancing product. This is the first time a PDE-5 inhibitor and a potential NO donor were identified in one molecule. A hydrolysis experiment showed that the new analogue was a prodrug of aildenafil and was therefore named **nitrosoprodenafil** (Venhuis et al., 2011) Fig 5(b)**.** Both PDE-5 inhibitors and nitrosamines cause vasodilatation by increasing levels of NO. To their coincidental use is warned against because it may cause a fatal drop in blood pressure. In addition, nitrosamines are known carcinogens. The findings indicate the dangerous level of advancement in medicinal chemistry by producers of unapproved drugs.

Tadalafil (Cialis®) was approved in 2003 by the FDA as the third phosphodiesterase type 5 enzyme (PDE-5) inhibitor to treat MED (Meuleman, 2003). Then, different tadalafil analogues have been found as adulterants in illegal products. Hasegawa et al. detected N**octyl-nortadalafil** Fig 6(a) together with cyclopentynafil **Fig** 6(b) in dietary supplement (Hasegawa et al., 2008). Both of them are the first compounds reported to be new tadalafil and sildenafil analogues.

O

N

N

Another sildenafil analogue was detected from a health supplement claimed for human MED, the structure of this new analogue was characterized as **dithiodesmethylcarbodenafil** containing 2 thiocarbonyl groups instead of 2 carbonyl groups, and 4-methyl substitution on the piperazine ring, rather than 4-ethyl substitution when

HN N

OH

O

N N

(b) (c)

HN

O

S O O

HN N

(a) (b)

Other new unapproved analogue of sildenafil was detected in capsules of a herbal dietary supplement promoted as a libido enhancing product. This is the first time a PDE-5 inhibitor and a potential NO donor were identified in one molecule. A hydrolysis experiment showed that the new analogue was a prodrug of aildenafil and was therefore named **nitrosoprodenafil** (Venhuis et al., 2011) Fig 5(b)**.** Both PDE-5 inhibitors and nitrosamines cause vasodilatation by increasing levels of NO. To their coincidental use is warned against because it may cause a fatal drop in blood pressure. In addition, nitrosamines are known carcinogens. The findings indicate the dangerous level of advancement in medicinal

Tadalafil (Cialis®) was approved in 2003 by the FDA as the third phosphodiesterase type 5 enzyme (PDE-5) inhibitor to treat MED (Meuleman, 2003). Then, different tadalafil analogues have been found as adulterants in illegal products. Hasegawa et al. detected N**octyl-nortadalafil** Fig 6(a) together with cyclopentynafil **Fig** 6(b) in dietary supplement (Hasegawa et al., 2008).

HN

N

chemistry by producers of unapproved drugs.

and sildenafil analogues.

 (a) Dithio-desmethylcarbodenafil; (b) **N**itroso-prodenafil.

O

Fig. 4. (a) Piperidenafil; (b) Piperazinonafil; (c) Isopiperazinonafil.

(a)

O

N N N

Fig. 6. (a) *N*-octylnor-tadalafil; (b) Cyclopentynafil.

In Taiwan, one of the dietary supplement which was claimed on the treatment of male erectile dysfunction was firstly screened in 2009 and Tadalafil and its doctored version was newly identified. Since it is having amino group instead of methyl in tadalafil it was named as **aminotadalafil** (Zou et al., 2006; Lin et al., 2009) Table-1. This compound has two asymmetric carbons, theoritically two pairs of enantiomers exist. The chromatographic separation of its stereoisomers was reported by using chiral LC-MS (Kurita et al., 2008). Using this method, *RR*-Aminotadalafil and *SR*-Aminotadalafil were detected in some health food. In addition, an interaction product of aminotadalafil was isolated from an illegal health food product. The structure of the interaction product was elucidated and unknown compound was characterized as condensation product of aminotadalafil and hydroxymethylfuraldehyde and is probably the result of a drug-excipient incompatibility (Häberli et al., 2010) Fig 7.

Fig. 7. Condensation product of aminotadalafil and hydroxymethylfuraldehyde.

Last flash development for the treatment of MED is discovering of **Zoraxel** (RX-10100) by Rexahn Pharmaceutical company (Albersen et al., 2010). Zoraxel is containing clavulanic acid that is centrally acting in the CNS and may be a more effective MED treatment for patients who are responsive or unresponsive to PDE-5 inhibitors. It is being developed as an orally administered, on-demand tablet to treat sexual dysfunction, and has extensive and well-established safety in humans. For the future, it is being expected, Zoraxel will be on to worldwide best-selling drug.

Chromatographic Separation and Identification of Sildenafil

\*Detected in our lab. for the first time by us.

Table 1. Some synthetic PDE-5 inhibitors as analogues sildenafil and tadalafil were isolated by us.

and Yohimbine Analogues Illegally Added in Herbal Supplements 57

In our central instrumental analysis laboratory, we also try to detect these commercially available supplements (which are sent by the Ministry of Food Agriculture and Livestock of Turkey before it grants a license for import to Turkey) whether they possess new or old sildenafil analogues by using high-performance liquid chromatography with diode array detection and mass spectrometry (HPLC-DAD–MS) and nuclear magnetic resonance (NMR) analyses, NMR is the only analytical technique which provides full structural information from novel compounds, acquisition of MS data or comparison of retention times may not be sufficient, in this case, LC-MS/NMR allowed to identify the adulterants without any need for references (Kesting et al., 2010) this means that these kind of analogues are not easy to detect by ordinary laboratory methods, we have identified the listed analogues are given in Table 1-2 up to date. With the aim of evaluating the potential risks of commercialized aphrodisiac products on consumer health, the aim of present work is to investigate simple HPLC-MS method and NMR data of synthetic and natural analogues of aphrodisiacs.

**Table-1** shows the determination of some sildenafil and tadalafil analogues and Dapoxetin HCl, in commercially available health supplements in Turkey. Their formulas, 1H-13C-NMR spectra, ESI(+) m/e values and their LC chromatograms are given in Table 1 and Fig 8, 9, respectively. Fig 9 shows the good separation of caffeine, some sildenafil analogues and

Fig. 8. HPLC chromatogram of *4.75* (Sildenafil), *4.80* (Vardenafil), *3.84* (Tadalafil), *3.62* (Aminotadalafil), *4.39* (Nor-carbodenafil) (*RT* : Retention times as min).

In our central instrumental analysis laboratory, we also try to detect these commercially available supplements (which are sent by the Ministry of Food Agriculture and Livestock of Turkey before it grants a license for import to Turkey) whether they possess new or old sildenafil analogues by using high-performance liquid chromatography with diode array detection and mass spectrometry (HPLC-DAD–MS) and nuclear magnetic resonance (NMR) analyses, NMR is the only analytical technique which provides full structural information from novel compounds, acquisition of MS data or comparison of retention times may not be sufficient, in this case, LC-MS/NMR allowed to identify the adulterants without any need for references (Kesting et al., 2010) this means that these kind of analogues are not easy to detect by ordinary laboratory methods, we have identified the listed analogues are given in Table 1-2 up to date. With the aim of evaluating the potential risks of commercialized aphrodisiac products on consumer health, the aim of present work is to investigate simple HPLC-MS method and NMR data of synthetic and

**Table-1** shows the determination of some sildenafil and tadalafil analogues and Dapoxetin HCl, in commercially available health supplements in Turkey. Their formulas, 1H-13C-NMR spectra, ESI(+) m/e values and their LC chromatograms are given in Table 1 and Fig 8, 9, respectively. Fig 9 shows the good separation of caffeine, some sildenafil analogues and

Sildenafil-T 2: Diode Array

Vardenafil-T 2: Diode Array

Tadalafil-T 2: Diode Array

Aminotadalafil-T 2: Diode Array

Nor-Carbodenafil-T 2: Diode Array

4.75

4.80

3.84

3.62

TIC 2.71e4

TIC 2.86e4

TIC 2.59e4

TIC 4.93e4

TIC 2.89e4

Fig. 8. HPLC chromatogram of *4.75* (Sildenafil), *4.80* (Vardenafil), *3.84* (Tadalafil), *3.62*

1.00 2.00 3.00 4.00 5.00 6.00 7.00 8.00 9.00 10.00 11.00 12.00 0 Time

(Aminotadalafil), *4.39* (Nor-carbodenafil) (*RT* : Retention times as min).

4.39

natural analogues of aphrodisiacs.

**Default file**

100 %

100 %

0

100 %

0

100 %

0

100 %

0

\*Detected in our lab. for the first time by us.

Table 1. Some synthetic PDE-5 inhibitors as analogues sildenafil and tadalafil were isolated by us.

Chromatographic Separation and Identification of Sildenafil

140.8 142.5 144.6 158.6 171.3 168.5

3.19

4.30

4.24

4.36

2.72

2.73 0.90 1.01 1.14 2.10 1.55 2.86 3.33 3.87

3.65 5.04 6.58 9.15 11.37 11.92 12.92

**Default file**

2.69

2.81 2.96 3.04 3.18

100

%

**Default file**

100

%

**Default file**

100

100

0

%

%

129.7 118.2

1.40 0.74

and Yohimbine Analogues Illegally Added in Herbal Supplements 59

L-ArgininHCl-T 1: Scan ES+

TIC 7.65e8

TIC 4.48e4

TIC 9.58e8

2.00 4.00 6.00 8.00 10.00 12.00 14.00 16.00 18.00 20.00 22.00 24.00 0 Time

L-ArgininHCl-T 285 (2.855) Cn (Cen,2, 80.00, Ht) 1: Scan ES+

175.5 2.57e7

212.8 222.5 238.5 207.8 177.7 187.8 189.5 197.4 236.9 251.5252.9 257.2

8.64 7.18

8.40 7.597.91

8.62

8.97

9.239.63 10.15 10.57 10.81 11.18

120 130 140 150 160 170 180 190 200 210 220 230 240 250 260 0 m/z

Natural 2: Diode Array

Natural 1: Scan ES+

7.18 4.28

4.584.80 5.13 5.51 6.49 5.92

176.4

Fig. 10. HPLC-MS ion chromatogram of *2.69* (L-Arginin HCl), *2.72* (Inosine), *3.19* (Icariin),

1.00 2.00 3.00 4.00 5.00 6.00 7.00 8.00 9.00 10.00 11.00 0 Time

*4.3* (Yohimbin), *7.18* (Imperatonin), *8.64* (Osthole). (*RT* : Retention times as min).

Fig. 9. HPLC chromotogram of *3.13* (Caffeine), *3.37* (Dimethylacetildenafil), *3.67* (Udenafil), *4.02* (Noracetildenafil), *4.43* (Dimethylsildenafil), *5.02* (Dimethylhomosildenafil), *6.43* (Hydroxythiohomosildenafil), *8.14* (Thiomethisosildenafil), *8.69* (Thiosildenafil), *9.24* (Dapoxetine HCl), *10.14* (Thiohomomethisosildenafil) (*RT* : Retention times as min).

Dapoxetin HCl. **Dapoxetine HCl**, is reported to be a selective serotonin reuptake inhibitor under investigation for the treatment of premature ejaculation (PE) (Li et al., 2009a). Since we have met this compound in one of the herbal drinks for adults, we have added it to our list. As it is well known, caffeine is very common for this kind of health supplements, it was added to chromatogram as well.

#### **1.1 Natural aphrodisiacs**

Formulas, 1H-13C-NMR spectra, ESI(+) m/e values and LC-MS ion chromatograms of natural aphrodisiacs are given in Table 2 and Fig. 10.

L- **Arginine** is not herb, but a nonessential amino acid, it is found naturally in foods such as meat, dairy, poultry and fish, it also may be synthesized in the laboratory, in spite of there is insufficient evidence to rate effectiveness for male fertility and female sexual problem, but it is possible to see the ingredients of many aphrodisiacs, recently. It is also available as oral Larginine supplements, which some product manufacturers market as a "natural Viagra" (Stanislavov & Nikolova, 2003).

1-9+Dap-Caff- 2: Diode Array

6.43 5.02 8.69

8.14

10.14

9.24

3.67

4.43 4.02

3.13

3.37

TIC 4.60e3

Fig. 9. HPLC chromotogram of *3.13* (Caffeine), *3.37* (Dimethylacetildenafil), *3.67* (Udenafil), *4.02* (Noracetildenafil), *4.43* (Dimethylsildenafil), *5.02* (Dimethylhomosildenafil), *6.43* (Hydroxythiohomosildenafil), *8.14* (Thiomethisosildenafil), *8.69* (Thiosildenafil), *9.24* (Dapoxetine HCl), *10.14* (Thiohomomethisosildenafil) (*RT* : Retention times as min).

1.00 2.00 3.00 4.00 5.00 6.00 7.00 8.00 9.00 10.00 11.00 12.00 0 Time

Dapoxetin HCl. **Dapoxetine HCl**, is reported to be a selective serotonin reuptake inhibitor under investigation for the treatment of premature ejaculation (PE) (Li et al., 2009a). Since we have met this compound in one of the herbal drinks for adults, we have added it to our list. As it is well known, caffeine is very common for this kind of health supplements, it was

Formulas, 1H-13C-NMR spectra, ESI(+) m/e values and LC-MS ion chromatograms of

L- **Arginine** is not herb, but a nonessential amino acid, it is found naturally in foods such as meat, dairy, poultry and fish, it also may be synthesized in the laboratory, in spite of there is insufficient evidence to rate effectiveness for male fertility and female sexual problem, but it is possible to see the ingredients of many aphrodisiacs, recently. It is also available as oral Larginine supplements, which some product manufacturers market as a "natural Viagra"

added to chromatogram as well.

(Stanislavov & Nikolova, 2003).

natural aphrodisiacs are given in Table 2 and Fig. 10.

**1.1 Natural aphrodisiacs**

**Default file**

100

%

Fig. 10. HPLC-MS ion chromatogram of *2.69* (L-Arginin HCl), *2.72* (Inosine), *3.19* (Icariin), *4.3* (Yohimbin), *7.18* (Imperatonin), *8.64* (Osthole). (*RT* : Retention times as min).

Chromatographic Separation and Identification of Sildenafil

currently recommended (Albersen et al., 2010).

and neurotoxicity in an Iranian male patient.

**2. Material and methods**

**2.1 Equipments**

supplements.

and Yohimbine Analogues Illegally Added in Herbal Supplements 61

**Yohimbine** is an alkaloid which are found naturally in *Pausinystalia yohimbe*, a classical aphrodisiac which has been recently revalued for its pro-sexual properties and extensively commercialized without control in some countries (Saini et al., 2010, Melnyk et al., 2011) . However, there is little evidence on its efficacy in the treatment of ED and it is, therefore, not

**Inosine** is a nucleoside consisting of ribose and hypoxanthine, most commonly found in supplements that claim to be "energy promoters", naturally it is found in brewer's yeast and organ meats, however it can be synthesized in laboratory (Deuster & Simmons, 2004). Their nerve-stimulating action probably enhances sexual functions or help to remedy sexual

**Icariin** is one of the primary active component of *Epimedium* extracts, which has been used to treat impotence and improve sexual function by acting as PDE5 inhibitor as nafil derivatives (Liu et al., 2005; Dell'Agli, 2008). Since the illegally using of synthetic nafil derivatives have been prohibited, recently icariin has become very attractive by the drug manufacturers and suppliers, because of its natural character. Nowadays, the percentage ratio of Icariin in the commercial *Epimedium* extracts has been increased over 80-90. However there is no sufficient information about its pharmacological profiles and safety. Further investigation needs to be done to examine any benefits that could occur from

**Osthole** and **imperatorin**, coumarin compounds have been reported to exhibit various biological activities (Liao et al., 2010). It was reported that, both of them were found to help relax the corpus cavernosa of the penis, which would potentially help with blood flow, in phenylephrine-precontracted endothelium-intact rabbit corpus cavernosum (Chen et al.,

Herbal medications are being progressively utilized all over the world and it is believed that herbal remedies are not hazards, however some adverse reactions have been increased. Tribulus terrestris is frequently used because of its aphrodisiac effect. But there is an article (Talasaz, 2010), which reports a case of *T-terrestris*-induced hepatotoxicity, nephrotoxicity

In this text, we also present a new sildenafil analogue was found to have been added illegally to a herbal drinks marketed for the enhancement of sexual function. This analogue has never been found. Therefore, it prompted us to elucidate its structure. The structure was determined as 5-[5-[[(3,5-dimethyl-1-piperazinyl]sulfonyl]-2-propoxyphenyl]-1,6-dihydro-1 methyl-3-propyl-7*H*-pyrazolo [4,3-d]pyrimidin-7-one. Owing to the inclusion of a methylene group in dimethylsildenafil, the detected compound was called **Dimethyhomosildenafil**. The sample was purified with column chromatography. The IR, HPLC-/MS (ESI+), and

Uncorrected melting points were measured on an Büchi B-540 capillary melting point apparatus. 1H (400 MHz) and 13C(100 MHz) NMR spectra were recorded employing a

completely assigned NMR data of **dimethylhomosildenafil** have been observed.

functions when these are reduced due to degeneration of nerve tissue.

2000; Chiou et al., 2001), however there is no information for human uses.

Table 2. Natural aphrodisiacs were isolated by us from some herbal supplements.

 Natural aphrodisiacs were isolated by us from some herbal supplements.

 However, there is little evidence on its efficacy in the treatment of ED and it is, therefore, not currently recommended (Albersen et al., 2010).

**Inosine** is a nucleoside consisting of ribose and hypoxanthine, most commonly found in supplements that claim to be "energy promoters", naturally it is found in brewer's yeast and organ meats, however it can be synthesized in laboratory (Deuster & Simmons, 2004). Their nerve-stimulating action probably enhances sexual functions or help to remedy sexual functions when these are reduced due to degeneration of nerve tissue.

**Icariin** is one of the primary active component of *Epimedium* extracts, which has been used to treat impotence and improve sexual function by acting as PDE5 inhibitor as nafil derivatives (Liu et al., 2005; Dell'Agli, 2008). Since the illegally using of synthetic nafil derivatives have been prohibited, recently icariin has become very attractive by the drug manufacturers and suppliers, because of its natural character. Nowadays, the percentage ratio of Icariin in the commercial *Epimedium* extracts has been increased over 80-90. However there is no sufficient information about its pharmacological profiles and safety. Further investigation needs to be done to examine any benefits that could occur from supplements.

**Osthole** and **imperatorin**, coumarin compounds have been reported to exhibit various biological activities (Liao et al., 2010). It was reported that, both of them were found to help relax the corpus cavernosa of the penis, which would potentially help with blood flow, in phenylephrine-precontracted endothelium-intact rabbit corpus cavernosum (Chen et al., 2000; Chiou et al., 2001), however there is no information for human uses.

Herbal medications are being progressively utilized all over the world and it is believed that herbal remedies are not hazards, however some adverse reactions have been increased. Tribulus terrestris is frequently used because of its aphrodisiac effect. But there is an article (Talasaz, 2010), which reports a case of *T-terrestris*-induced hepatotoxicity, nephrotoxicity and neurotoxicity in an Iranian male patient.

In this text, we also present a new sildenafil analogue was found to have been added illegally to a herbal drinks marketed for the enhancement of sexual function. This analogue has never been found. Therefore, it prompted us to elucidate its structure. The structure was determined as 5-[5-[[(3,5-dimethyl-1-piperazinyl]sulfonyl]-2-propoxyphenyl]-1,6-dihydro-1 methyl-3-propyl-7*H*-pyrazolo [4,3-d]pyrimidin-7-one. Owing to the inclusion of a methylene group in dimethylsildenafil, the detected compound was called **Dimethyhomosildenafil**. The sample was purified with column chromatography. The IR, HPLC-/MS (ESI+), and completely assigned NMR data of **dimethylhomosildenafil** have been observed.

#### **2. Material and methods**

#### **2.1 Equipments**

Uncorrected melting points were measured on an Büchi B-540 capillary melting point apparatus. 1H (400 MHz) and 13C(100 MHz) NMR spectra were recorded employing a

Chromatographic Separation and Identification of Sildenafil

Fig. 11. IR spectrum of Dimethylhomosildenafil.

toxicity profile.

and Yohimbine Analogues Illegally Added in Herbal Supplements 63

methisosildenafil (Reepmeyer & Avignon, 2009). The IR spectrum of dimethylhomosildenafil is given in Fig 11. This compound also must be put on the inspection list for illegal health-related substances because of the unknown safety and

VARIAN MERCURY 400 MHz FT spectrometer, with CDCl3 as solvent. Chemical shifts (δ) are in ppm relative to TMS. The LC/MS were taken on a Waters Micromass ZQ connected with Waters Alliance HPLC, using ESI(+) method, with C-18 column. Elemental analyses were performed by Leco CHNS-932. The infrared spectrum was recorded in the 600-3600 cm-1 range using a Jasco FT-IR-420 spectrometer and KBr pellets.


**1. Introduction**

154 Rheumatoid Arthritis – Treatment

Zwerina J, Hayer S, Tohidast-Akrad M, Bergmeister H, Redlich K, Feige U, et al. Single and

erosion, and cartilage destruction. Arthritis Rheum. 2004;50(1):277-90.

combined inhibition of tumor necrosis factor, interleukin-1, and RANKL pathways in tumor necrosis factor-induced arthritis: effects on synovial inflammation, bone

> Enkorten is a new potential drug for the treatment of rheumatoid arthritis**,** with an immunomodulatory and anti-inflammatory effect. It is a combination of two peptide components of endogenous origin: methionine-enkephalin of 5 mg and tridecactide of 1 mg (Picture 1 and 2). According to the chemical structures, these components correspond to amino acid sequences of the neuropeptide precursor proopiomelanocortin.

Picture 1. Chemical characteristics of tridecactide

Picture 2. Chemical characteristics of met-enkephalin

The holder of effects of the tested combination is alpha 1-13 corticotrophin–tridecactide (ACTH 1-13), composed of 13 amino acids, identical to the amino acid sequence of the alpha-melanostimulating hormone (-MSH), but without the acetyl and amide ending. Some studies suggest that -MSH is synthesized in the pituitary, other parts of the CNS, the placenta, and some endocrine organs – especially the adrenal, the skin and the gastrointestinal tract (Star et al, 1995; Catania et al, 2000). The presence of this neuropeptide in the mucous membrane of the gastrointestinal tract and keratinocytes of the skin suggests

<sup>\*</sup> Jasna Kusturica, Maida Todić-Rakanović, Mirjana Mijanović, Fahir Bečić, Asija Začiragić,

Selma Škrbo, Lejla Burnazović-Ristić, Aida Kulo and Svjetlana Loga-Zec

*Medical Faculty University of Sarajevo Institute of Pharmacology, Clinical Pharmacology and Toxicology, Bosnia and Herzegovina*

Enkorten – A Potential Drug for the Treatment of Rheumatoid Arthritis 157

inflammatory response of T cells through the activation of regulatory physiological

Until recently it was thought that the immune/inflammatory reaction develops exclusively in the periphery. Today it is known that the nervous, endocrine and immune systems are interconnected, and that they communicate by signals that are transmitted through neuropeptides. In this way, central neurogenic influences can either enhance or modulate the peripheral response, so by analogy the treatment of inflammatory reactions could be improved by the influence of the proinflammatory signals in the central nervous system. Thus, Catania et al. (1999) suggest that one of these strategies could be based on the action of

In the past several years, particular attention has been given to the anti-inflammatory effects of this peptide, as a result of the introduction of -MSH and related peptides in clinical practice for the treatment of inflammatory diseases (Getting, 2006). Undoubtedly, -MSH, administered systemically or locally, expresses strong anti-inflammatory effects. Antiinflammatory effects of this peptide are mediated by direct effects on cells of the immune system, as well as by indirect effects that are achieved by changing the function of non-

Three basic mechanisms of -MSH anti-inflammatory action can be distinguished: direct effects through melanocortin receptors on cells in the periphery (monocytes/macrophages and neutrophils); effects on glial cells; and descendent anti-inflammatory effects through melanocortin receptors on neurons (Lipton et al, 2000). These mechanisms of action correspond to the concept of neuroimmunomodulation in which nerves, endocrine and immune systems are independent networks that mutually communicate through soluble

Multiple anti-inflammatory effects that -MSH has have been identified in various *in vitro* systems of cell cultures. Studies of anti-inflammatory effects of -MSH on cellular level were primarily concentrated on the issue whether and to what degree -MSH suppresses the production of pro-inflammatory cytokines. In this respect, several studies have shown that treatment with -MSH results in significant downregulation of pro-inflammatory cytokines such as IL-1, IL-2, IL-4, IL-6, IL-13 and TNF-, as well as chemokines such as IL8, Gro and interferon (IFN) (Brzoska et al, 2008). Furthermore, chemotaxis induced by IL8 in the cells of human neutrophils and monocytes has been suppressed under the influence of -MSH (Manna et al, 2006), which indicates that the function of these types of phagocytic cells during the inflammatory response is blocked by the peptide via multiple effector pathways. Contrary to the inhibitory effects of -MSH on the production and activity of proinflammatory mediators, it has been found that -MSH induces IL-10, a cytokine with a potent anti-inflammatory activity. Specifically, in monocytes of human peripheral blood and cultivated human monocytes, -MSH increases the production and expression of IL-10 (Bhardway et al, 1996). Since IL-10 reduces the production of pro-inflammatory cytokines in

resident immune cells of peripheral tissues (Luger and Brzoska, 2007).

macrophages, its upregulation could result in anti-inflammatory influences.

Studies on different cell types of human skin, including pigment cells, fibroblastic cells and dermal microvascular endothelial cells, as well as rat mast cells, showed that -MSH is able to suppress the expression of several intracellular adhesion molecules (interstitial adhesion molecule-1 - ICAM-1 and P-selectin) induced by proinflammatory stimuli such as IFN-, LPS and TNF-Böhm et al, 2005; Hill et al, 2006; Sarkar et al, 2003). Finally, the inhibition of adhesion and transmigration of inflammatory cells may contribute to the anti-

mediators such as cytokines and neuropeptides.

mechanisms (Štambuk et al, 1997).

that it is a component of natural immunity (Catania et al, 2000a). The supposed connection between alpha-MSH and immune processes is confirmed by the clinical data on changes of its plasma level under certain pathological conditions (Lipton and Catania, 1997; Lipton et al, 2000). An increased concentration of -MSH has been recorded in patients with an acute myocardial infarction who received thrombolytic therapy, probably as an endogenous antiinflammatory response. This peptide has been identified in the synovial fluid of patients with rheumatoid arthritis (Star et al, 1995) and it is believed that the inflammation and joint destruction in this disease occur due to imbalances of production of proinflammatory and anti-inflammatory cytokines such as -MSH. A higher concentration of -MSH has been found in patients suffering from arthropathy accompanied by heavy inflammation (Lipton et al, 2000). Since -MSH is a natural modulator peptide, its concentration is subject to intensive control to avoid excessive inhibition of the immune response. It is assumed that the increased concentration of peptides leads to downregulation of receptors in cells, which could explain the biphasic curve of the anti-inflammatory effects -MSH (Catania et al, 2000).

The other component contained in the drug is met-enkephalin, which belongs to a group of endogenous opioid peptides, widely distributed in the CNS. All opioid peptides contain one or more copies of the simplest opioid peptides: methionine-enkefalin (tyrosine-glycineglycine-phenylalanine-methionine) and leucine-enkephalin (tyrosine-glycine-glycinephenylalanine-leucine) (Janković, 1991). They participate in, and modulate neurotransmission. Opioid peptides also act as growth factors that control proliferation and differentiation of cells and may participate in the process of wound healing, tissue regeneration and immune response. Enkephalins show antidepressant, antianxiety and anticonvulsive effects. They are produced in the adrenal and the hypothalamus. Metenkephalin belongs to a group of enkephalins that are involved in transmitting signals in the nervous system, intestines, endocrine glands, bronchial tree, lymphatic and other tissues (Štambuk et al, 1996). The highest concentrations are present in the core of adrenal, hypothalamus, autonomic ganglia and the gastrointestinal tract (Janković, 1991).

Met-enkephalin as an inhibitory peptide modulates the proliferation and migration of cells, as well as the organization of tissues during development, homeostatic cell renewal, wound healing, angiogenesis, and in malignant neoplasms (Zagon et al, 2000). Physiologically, enkephalins have a signalling function in the nervous system, intestine, endocrine glands, bronchial tree, cardiovascular system, lymphatic and other tissues, and participate in the regulation of haematopoiesis (Štambuk et al, 1996). It is assumed that they participate in maintaining homeostasis of the organism exposed to stress (Mulabegović et Rakanović-Todić, 2008).