## 1 Introduction  
Contemporary antenna structures are designed to meet the demands pertaining to different application areas, including wireless communications [1] (along with the emerging 5G technology [2, 3]), satellite communications [4], medical imaging [5], or internet of things (IoT) [6]. In many cases, various functionalities have to be implemented including multi-band operation [7], circular polarization [8], band notches [9], MIMO operation [10], or pattern diversity [11]. At the same time, reduction of the physical size of radiators is of concern for a growing number of situations (e.g., wearable [12] and implantable devices [13]). In pursuit of fulfilling the stringent specifications, increasingly complex antenna topologies are being developed. Their evaluation requires full-wave electromagnetic (EM) analysis because other means are either unavailable or (as equivalent network models) can only be used for yielding rough initial designs that need further tuning [14]. Nowadays, EM-driven design is ubiquitous, but its downside is high computational cost. This becomes a practical issue whenever multiple simulations are necessary. A typical task is parameter adjustment [15; 16]. It is imperative for most antenna structures but generates considerable CPU expenditures, particularly when reasonable initial design is not available, or the global search needs to be involved [17]. Expediting design optimization can be realized at the algorithmic level. Representative examples are the gradient-based procedures where the speedup can be achieved by lowering the cost of antenna response Jacobian estimation using adjoint sensitivities [18; 19], or by suppressing finite-differentiation sensitivity updates based on various criteria, e.g., monitoring design relocation [20] or gradient changes between the algorithm iterations [21]. The employment of fast surrogate models provides viable alternatives. These can be realized in the context of local [22] or global optimization [23]. For the former, the surrogates are typically constructed along the optimization path, and can be data-driven (polynomial regression [24], kriging [25], neural networks [26], polynomial chaos expansion [27]) or physics-based (space mapping [28; 29], response correction [30], feature-based optimization [31]), normally obtained from underlying low-fidelity models (in the case of antennas most often being coarse-mesh EM simulations [32]). Global surrogate-assisted optimization often involves machine learning techniques [33] but also all kinds of approximation surrogates [34]. In this context, sequential sampling methods are typically employed [35], where iterative allocation of the infill points can be aimed at improving the predictive power of the model, searching for the global optimum or the combination thereof [36]. Utilization of the stand-alone surrogate models as replacements of CPU-heavy EM simulations has been fostered in the literature as a way of accelerating various simulation-driven design procedures, including design closure [37; 38; 39]. Data-driven surrogates seem to be especially attractive for these purposes due to their versatility and widespread availability (e.g., [40; 41]). Well established and popular approximation techniques include polynomial regression [36], artificial neural networks [42], radial basis functions (RBF) [37], kriging interpolation [36], Gaussian process regression (GPR) [38], or support vector regression (SVR) [43]. Unfortunately, conventional methods are severely limited when applied to antenna problems. The primary challenges include nonlinearity of antenna characteristics, the curse of dimensionality affecting the modeling process due to a typically large number of geometry/material parameters that need to be handled, as well as the need for covering wide ranges of operating conditions and designable parameters. The latter is critical from the point of view of practical usefulness of the surrogate for design purposes.