## 3 SPSA Learning-Based Three-Dimensional Fuzzy Modeling
### SPSA Learning for Resulting Components of Three-Dimensional Fuzzy Rule
#### 3.3.2 Parameter Learning Using SPSA Algorithm  
SPSA [24] is an optimization algorithm based on stochastic approximation. In contrast with the gradient descent method, SPSA does not necessitate the direct measurement of the gradient of the target function. Instead, it just requires two measurements of the target function at every iteration to approximate the gradient, regardless of the complexity of the optimization problem. The two measurements are obtained by randomly and simultaneously altering all the variables in the problem in a correct manner, which is known as \"simultaneous perturbation\". SPSA is particularly suitable for nonlinear optimization problems with high-dimensional parameter space, which are difficult to calculate the gradient or have no analytical expression.  
In this study, in each iteration, the approximation gradients of all parameters of Fourier space base functions are computed by two randomly disturbed measurements of target functions. Let \(\hat{\vartheta}_{k}=[a_{1,k},b_{1,k},c_{1,k},d_{1,k},a_{2,k},b_{2,k},c_{2,k},d_{ 2,k},d_{2,k},\ldots,a_{Z^{\prime},k},b_{Z^{\prime},k},c_{Z^{\prime},k},d_{Z^{ \prime},k}]\in R^{4Z^{\prime}}\) be all the parameters at the \(k\)th iteration, where \(a_{i,k},b_{i,k},c_{i,k}\) and \(d_{i,k}\) are parameters of the \(i\)th space base function at the \(k\)th iteration. The target function is described as the L2 norm of the difference between the spatial output \(w(\bar{F},t)\) and the predicted spatial output \(\hat{w}(\bar{F},t)\) as follows:  
\[J(\hat{\vartheta})=\frac{1}{2}\sum_{i=1}^{n}(w_{i}(\bar{F},t)-\bar{w}_{i}( \bar{F},t))^{2}. \tag{12}\]  
The algorithm of SPSA is described as follows.  
**Step 1: Initialization and selection of coefficients**  
Initialize the timer number as \(k=0\). Select an initial estimate \(\hat{\vartheta}_{0}\) and choose non-negative constants \(b\), \(d\), \(A\), \(\alpha\), and \(\beta\) for the SPSA algorithm. Then, calculate sequences \(b_{k}\) and \(d_{k}\) as in (13) and (14):  
\[b_{k}=b/(A+k+1)^{a}. \tag{13}\]  
\[d_{k}=d/(k+1)^{\beta}. \tag{14}\]  
**Step 2: Simultaneous perturbation vector generation**  
Create a random perturbation vector \(\Delta_{k}\) in a \(4Z^{\prime}\)-dimensional space, where each component of \(\Delta_{k}\) is individually produced from a zero-mean probability distribution that satisfies the specified constraints in [24]. Here, we use Bernoulli \(\pm 1\) distribution, which means that there is an equal probability of 0.5 for each event of \(\pm 1\).  
**Step 3: Evaluation of the cost function**  
Calculate two values of the cost function \(J(\hat{\vartheta})\) by perturbing the current \(\hat{\vartheta}_{k}\) simultaneously: \(J\left(\hat{\vartheta}_{k}+d_{k}\Delta_{k}\right)\) and \(J\left(\hat{\vartheta}_{k}-d_{k}\Delta_{k}\right)\), where \(d_{k}\) and \(\Delta_{k}\) are obtained from Step 1 and Step 2, respectively.  
**Step 4: Approximation of the gradient**  
Compute the simultaneous perturbation estimate for the gradient \(\mathbf{g}(\hat{\vartheta}_{k})\):  
\[\hat{\mathbf{g}}_{k}\big{(}\hat{\vartheta}_{k}\big{)}=\frac{\int\big{(}\hat{ \vartheta}k+d_{k}\Delta_{k}\big{)}-\int\big{(}\hat{\vartheta}_{k}-d_{k}\Delta_ {k}\big{)}}{2d_{k}}\left[\begin{array}{c}\Delta_{k1}\\ \Delta_{k2}\\ \vdots\\ \Delta_{kp}\end{array}\right]. \tag{15}\]  
where \(\Delta_{ki}\) represents the \(i\)th element of the \(\Delta_{k}\) vector.  
**Step 5: Revision of the \(\vartheta\) estimation**  
\[\hat{\vartheta}_{k+1}=\hat{\vartheta}_{k}-b_{k}\hat{\mathbf{g}}_{k}\big{(} \hat{\vartheta}_{k}\big{)}. \tag{16}\]  
Using (16) to modify \(\hat{\vartheta}_{k}\) to a fresh value \(\hat{\vartheta}_{k+1}\).  
**Step 6: Iterating or terminating**  
Go back to Step 2 and substitute \(k\) with \(k+1\). Terminate the algorithm if there is minimal variation in several consecutive iterations or if the maximum number of iterations has been reached.