Fractals in Neural Networks: Introducing a New Class of Activation Functions

1. Introduction

Feed-forward artificial neural networks (ANNs) approximate a mapping $x\mapsto y$ by stacking linear layers interleaved with element-wise activation functions. While the weights learn a task-specific projection, the activations $\varphi (·)$ inject the non-linearity that lets the model fit highly curved decision boundaries. Classical choices range from smooth, saturating [1,2,3] (sigmoid, tanh) to piece-wise linear [4,5] (ReLU) functions. Deeper networks with carefully chosen $\varphi$ achieve state-of-the-art results in vision, language, and scientific discovery. This success is closely related to the expressive power of the activations. This expressiveness depends on both the choice of the activation function and on the shape of $\varphi$ [6,7].

Recent evidence  [8] indicates that the boundary separating trainable and divergent hyper-parameter configurations is itself a fractal object. The update map iterated during training behaves much like the complex quadratic map that generates Mandelbrot and Julia fractals—tiny changes in learning rate or initialisation can flip the outcome from convergence to blow-up.

We took this emergence of fractality within neural networks as inspiration to ask: What happens if we build fractality into neural networks from the start?

Standard sigmoid-based shallow nets are limited in expressivity; greater depth and activations such as ReLU improve expressivity [6,7]. Here we introduce fractal activation functions to push expressivity beyond ReLU. Our designs draw on the Blancmange (Takagi) curve [9,10] and Weierstrass–Mandelbrot functions [11,12,13,14], blended with familiar components such as ReLU or tanh.

However, at this point we need to differentiate our fractal activation functions from other approaches engineered fractality in neural networks, as presented by [15,16,17,18,19]. These approaches are based on different concepts: while we introduce fractal activation functions for neural networks, they employ a fractal hierarchy of neurons and fractal features as a preprocessing step, respectively. Also similar to these approaches, there is evidence that fractal interpolation can improve a neural network’s performance on specific datasets [20,21]. Fractal techniques can also be used to improve convolutional neural network architectures and their hyperparameters [22,23].

This work addresses two research questions:

• Can activation functions derived from fractal series be used in neural networks like standard activations, and how do they affect performance?
• Do such fractal activations increase the expressive capacity of neural networks compared with existing choices?

Contributions.

(i)

We present a general recipe for converting Weierstrass- and Blancmange-type series into computationally stable activation functions (Section 2.3) and demonstrate their usability on ten public classification datasets.

(ii)

We quantify their expressivity through trajectory-length diagnostics, revealing super-ReLU growth and distinctive oscillatory fingerprints (Section 4).

Taken together, our results suggest that engineered fractality is a promising new paradigm —additionally to depth, width, and regularization—for improving the expressiveness and performance of neural networks.

2. Methodology

This section states the computational tools used in our study: feed-forward neural networks, five gradient-based optimizers, three standard activation functions and fractal functions.

2.1. Neural Networks

A feed-forward neural network realises a composite function $f(\mathbf{x};\theta ):{\mathbb{R}}^d\to {\mathbb{R}}^k$ by stacking L affine layers and point-wise non-linearities. For layer $\ell \in \{1,...,L\}$ with weight matrix $W^{(\ell )}\in {\mathbb{R}}^{m_{\ell }\times n_{\ell }}$ and bias vector ${\mathbf{b}}^{(\ell )}\in {\mathbb{R}}^{n_{\ell }}$,

$${\mathbf{z}}^{(\ell )}=W^{(\ell )}{\mathbf{a}}^{(\ell -1)}+{\mathbf{b}}^{(\ell )},{\mathbf{a}}^{(\ell )}=\varphi \left({\mathbf{z}}^{(\ell )}\right),$$

where ${\mathbf{a}}^{\left(0\right)}=\mathbf{x}$ and $\theta ={\{W^{(\ell )},{\mathbf{b}}^{(\ell )}\}}_{\ell =1}^L$. Gradients ${\nabla }_{\theta }\mathcal{L}$ of a differentiable loss $\mathcal{L}$ are computed by back-propagation [1] and supplied to an optimizer.

Parameter updates use one of the following first-order methods:

• Stochastic Gradient Descent (SGD) [24]
• RMSProp [25]
• Adam [26]
• Adagrad [27]
• Adadelta [28]

All optimizers adapt the raw gradient to improve convergence speed or stability; their update rules are given in the cited sources.

In this research we use three standard activation functions, [4,5]:

$$\begin{array}{cc}\hfill \mathrm{ReLU}\text{:}& \varphi \left(z\right)=max(0,z)[4,5],\hfill \\ \hfill \mathrm{Sigmoid}\text{:}& \varphi \left(z\right)={\displaystyle \frac{1}{1+e^{-z}}}[1,2],\hfill \\ \hfill \mathrm{Tan}\mathrm{h}\text{:}& \varphi \left(z\right)=tanh\left(z\right)={\displaystyle \frac{e^z-e^{-z}}{e^z+e^{-z}}}\left[3\right].\hfill \end{array}$$

2.2. Fractal Functions

Conventional activation functions (ReLU, tanh, sigmoid) create either a fixed number of linear pieces (ReLU) or a single sigmoidal bend (tanh/sigmoid). A fractal activation instead super-imposes infinitely many self-similar ripples of shrinking amplitude and growing frequency. Formally, one starts with a periodic function $\phi \left(x\right)$ and builds a weighted sum

$$\mathcal{F}\left(x\right)=\sum _{n=0}^{\infty }{a}^n\phi \left(b^{n}x\right),0<a<1,b>1.$$

(1)

The geometric weights $a^n$ guarantee uniform convergence for uniformly continuous $\phi$, meaning that $\mathcal{F}$ is continuous, while the ever finer oscillations imply that $\mathcal{F}$ is nowhere differentiable in the classical sense for a lot of choices of $\phi \left(x\right)$ — the hallmark of fractal graphs. Nowhere differentiability means that the function does not admit a finite derivative at any point in its domain, despite being continuous everywhere. This property has been rigorously studied for Weierstrass-type functions; see, for example, [14]. In practice we truncate the sum at $N\approx 30-100$, the higher harmonics are already below machine precision.

Why use such curves in a network?

• They inject high-frequency detail even in shallow layers, increasing functional richness without extra parameters.
• They retain bounded outputs and controlled Lipschitz constants when $a<1/b$ and $\phi$ is Lipschitz continuous, [14].
• Most importantly, for gradient-based optimization, automatic differentiation remains well defined through finite truncation used in code.

The Blancmange or Takagi function $\mathcal{B}$,9,10] replaces the base wave $\phi$ in (1) with the saw-tooth function, which calculates the absolute distance of x to the nearest integer

$$s\left(x\right)=dist\left(x,\mathbb{Z}\right)=\underset{k\in \mathbb{Z}}{min}|x-k|.$$

We further modulate its slope by a smoothing gate $tanh\left(\beta x\right)$ to avoid excessively flat regions near zero:

$$\mathcal{B}\left(x\right)=\sum _{n=0}^{\infty }{\displaystyle \frac{tanh\left(\beta 2^{n}x\right)}{2^n}}s\left(2^{n}x\right),\beta =1.$$

(2)

$\mathcal{B}$ is continuous, 1-periodic, bounded between 0 and $1/2$, odd around  $x=0$, and exhibits a dense set of corner-points that provide sharp gradient signals during training.

The Weierstrass function $W_{\alpha ,\beta }\left(x\right)={\sum }_{n=0}^{\infty }{\alpha }^{n}cos\left({\beta }^n\pi x\right)$ is the prototypical example of a continuous, yet nowhere differentiable function. We adapt it in two ways:

a)

Exponential envelope. Multiplying by $e^{-\lambda \left|x\right|}$localises the oscillations and keeps activations bounded for large $\left|x\right|$.

b)

Non-linear carriers. Replacing the cosine by $x{sin}^2(·)$, $tanh\left(x\right)+sin(·)$ or $ReLU\left(x\right)+sin(·)$ injects asymmetric slopes that interact well with SGD-like optimisers.

With $0<\alpha <1$, $\beta >1$ and $\lambda >0$

$${\mathcal{W}}_{\mathrm{env}}\left(x\right)=\sum _{n=0}^{\infty }{(-1)}^n{\alpha }^n\phi \left({\beta }^n\pi x\right)e^{-\lambda \left|x\right|}$$

(3)

where $\phi (·)$ stands for one of

$$cos(·),x{sin}^2(·),tanh\left(x\right)+sin(·),ReLU\left(x\right)+sin(·).$$

For all choices above ${\mathcal{W}}_{\mathrm{env}}$ is Lipschitz,1 odd, and oscillatory with an effective period that shrinks from 2 at small $\left|x\right|$ to almost 0 as $\left|x\right|\to \infty$.

2.3. Fractal Activation Functions

Each activation below is the finite partial sum of an otherwise infinite fractal series. The full series would be continuous yet nowhere differentiable; its classical derivative vanishes almost everywhere, making back-propagation inoperable. We therefore truncate the sum of the series at a fixed number of terms N, preserving the self-similar oscillations on the scales most relevant in practice while keeping the functions fully differentiable.2

This results in the following activation functions:

2.3.1. Modulated Blancmange Curve blancmange

$$\begin{array}{cc}\hfill f\left(x\right)& =\sum _{n=0}^{N-1}tanh\left(a{2}^{n}x\right){\displaystyle \frac{|x{2}^{n}mod2-a\sqrt{\left|x\right|}|}{2^n}}\hfill \end{array}$$

(4)

$$\begin{array}{cc}\hfill & a=0.75,N=30.\hfill \end{array}$$

(5)

2.3.2. Decaying Cosine Function decaying_cosine

$$\begin{array}{cc}\hfill g\left(x\right)& =\sum _{n=0}^{N-1}\zeta \left(0.05tanh\left(\pi x\right)+c^{n}cos\left(d^n\pi x\right)e^{-\left|x\right|/2}sign\left(x\right)\right)\hfill \end{array}$$

(6)

$$\begin{array}{cc}\hfill & c=0.5,d=2,\zeta =0.2666,N=75.\hfill \end{array}$$

(7)

2.3.3. Modified Weierstrass–Tanh modified_weierstrass_tanh

$$\begin{array}{cc}\hfill j\left(x\right)& =\sum _{n=0}^{N-1}{(-1)}^n{a}^{n}cos\left(b^n\pi x\right)e^{-0.75\left|x\right|}+tanh\left(x\right)\hfill \end{array}$$

(8)

$$\begin{array}{cc}\hfill & a=0.5,b=1.5,N=100.\hfill \end{array}$$

(9)

2.3.4. Modified Weierstrass–ReLU modified_weierstrass_ReLU

$$\begin{array}{cc}\hfill k\left(x\right)& =\sum _{n=0}^{N-1}{(-1)}^n{a}^{n}cos\left(b^n\pi x\right)e^{-0.75\left|x\right|}+ReLU\left(x\right)\hfill \end{array}$$

(10)

$$\begin{array}{cc}\hfill & a=0.5,b=3,N=100.\hfill \end{array}$$

(11)

2.3.5. Classical Weierstrass Function weierstrass

$$\begin{array}{cc}\hfill w\left(x\right)& =\sum _{k=1}^N{2}^{-k\gamma }sin\left(2\pi {\lambda }^{k}x\right)\hfill \end{array}$$

(12)

$$\begin{array}{cc}\hfill & \gamma =0.5,\lambda =2,N=100.\hfill \end{array}$$

(13)

This “plain” Weierstrass curve serves as a sanity-check baseline against its modulated counterparts.

2.3.6. Weierstrass–Mandelbrot Variants

$$\begin{array}{cc}\hfill m\left(x\right)& =\sum _{k=1}^N{2}^{-k\gamma }x{sin}^2\left(2\pi {\lambda }^{k}x\right),\hfill \\ \hfill & \mathtt{weierstrass}\_\mathtt{mandelbrot}\_\mathtt{xsinsquared}\hfill \\ \hfill n\left(x\right)& =\sum _{k=1}^N{2}^{-k\gamma }\left(x+sin\left(2\pi {\lambda }^{k}x\right)\right),\hfill \\ \hfill & \mathtt{weierstrass}\_\mathtt{mandelbrot}\_\mathtt{xpsin}\hfill \\ \hfill o\left(x\right)& =\sum _{k=1}^N{2}^{-k\gamma }\left(ReLU\left(x\right)+sin\left(2\pi {\lambda }^{k}x\right)\right),\hfill \\ \hfill & \mathtt{weierstrass}\_\mathtt{mandelbrot}\_\mathtt{relupsin}\hfill \\ \hfill p\left(x\right)& =\sum _{k=1}^N{2}^{-k\gamma }\left(tanh\left(x\right)+sin\left(2\pi {\lambda }^{k}x\right)\right),\hfill \\ \hfill & \mathtt{weierstrass}\_\mathtt{mandelbrot}\_\mathtt{tanhpsin}\hfill \end{array}$$

with $\gamma =0.5,\lambda =2,N=100$.

Naming convention. The above names are used throughout the code listings, tables, and figures (e.g. Table 3 lists decaying_cosine, modified_weierstrass_tanh, ...).

Figure 1. Mixed sets of activation functions over the interval $[-2,2]$. Each curve is drawn with the Tailwind-inspired five-colour palette and distinct line styles for clarity. Classical (ReLU, sigmoid, tanh) and fractal variants are deliberately interleaved to highlight differences in slope saturation, oscillations, and piece-wise linear segments.

Figure 1. Mixed sets of activation functions over the interval $[-2,2]$. Each curve is drawn with the Tailwind-inspired five-colour palette and distinct line styles for clarity. Classical (ReLU, sigmoid, tanh) and fractal variants are deliberately interleaved to highlight differences in slope saturation, oscillations, and piece-wise linear segments.

3. Classification Experiments

The primary aim of this research is to demonstrate the practical viability and effectiveness of our proposed methodology, thereby serving as a proof of concept for the underlying ideas. We test our previously discussed ideas on a diverse collection of publicly accessible datasets, all of them classification tasks.

Furthermore, we partition each dataset using an 80/20 split for training and testing, respectively. This division allows for extensive model training while also providing a robust evaluation of model performance on unseen data.

The effectiveness of our algorithms is shown using accuracy as the primary metric. Furthermore, we do not rely solely on the accuracy of a single run; instead, we perform 40 different train/test splits with an 80/20 ratio for each dataset, training the employed neural networks separately for each split. Thus, our experiments do not yield just a single prediction result or accuracy value, but rather a distribution that reflects how well the different activation functions performed for a given problem. For each dataset, we chose a fixed neural network architecture, for which we varied the optimizer and the activation function, resulting in a total of 60 different configurations/results per dataset (5 different optimizers × 12 tested activation functions). We also kept the initialization of the weights of the network the same for all different train/test splits. The architecture was determined by first training a neural network with RMSProp and a ReLU activation function to obtain reasonable baseline results. Subsequently, we systematically interchanged the activation functions and optimizers for the actual experiments.

The following details the employed datasets and the outcomes of all runs; Section 3.2 (Table 1, Table 2, Table 3, Table 4, Table 5, Table 6, Table 7, Table 8, Table 9 and Table 10)

3.1. Datasets

We employed the following classification task data sets, and neural network architectures to verify our ideas. For every dataset we fix the neural–network hyper-parameters once and keep them identical across all optimizers and activation functions.

Climate-Model Simulation Crashes Dataset

(OpenML: climate-model-simulation-crashes, ID: 1467):

540 simulations characterised by 18 numerical configuration variables (temperature, pressure, humidity, …). Target = crash (1) vs. no crash (0). Binary failure prediction of climate runs. Network: one hidden layer with 32 neurons and a sigmoid output; 30 epochs, batch size 32.

Diabetes Dataset

(OpenML: diabetes, ID: 37):

768 instances, 8 numeric predictors taken from diagnostic measurements (glucose, blood pressure, BMI, …). Target indicates the presence (1) or absence (0) of Type-2 diabetes in female patients of Pima Indian heritage. Binary classification. Network: one hidden layer with 64 neurons and a 1-unit sigmoid head; 30 epochs, batch size 32.

Glass Identification Dataset

(OpenML: glass, ID: 41):

214 instances described by 9 continuous chemical/physical attributes (e.g. refractive index, Na, Mg, Al, … % weight). The goal is to classify each sample into one of 6 glass–type categories (building windows float, building windows non–float, vehicle windows, containers, tableware, headlamps). A 6-class classification dataset. Network: two hidden layers with 64 and 128 neurons, followed by a 6-way soft-max output; trained for 30 epochs with mini-batches of 32.

Ionosphere Dataset

(OpenML: ionosphere, ID: 59):

351 radar returns (instances) with 34 real-valued attributes that summarize the complex temporal signals emitted in the ionosphere experiment at Lincoln Laboratory. The task is to predict good vs. bad radar returns, making this a binary classification problem. Network: one hidden layer with 128 neurons plus a single sigmoid output; 30 epochs, batch size 32.

Iris Dataset

(OpenML: iris, ID: 61):

150 flower samples with 4 continuous measurements (sepal length/width, petal length/width). Target species: setosa, versicolor, virginica. Classic 3-class benchmark dataset. Network: two hidden layers with 32 and 64 neurons; 3-way soft-max output; 30 epochs, batch size 16.

Seeds Dataset

(OpenML: seeds, ID: 1499):

210 grain kernels measured via 7 geometric descriptors (area, perimeter, compactness, length, width, asymmetry, groove length). The aim is to classify them into 3 wheat varieties (Kama, Rosa, Canadian). A 3-class problem with only numeric features. Network: two hidden layers with 64 and 128 neurons; 3-way soft-max output; 30 epochs, batch size 32.

Tic-Tac-Toe Endgame Dataset

(OpenML: tic-tac-toe, ID: 50):

958 board positions (instances), each encoded by 9 categorical attributes (cell contents {x,o,b}). Task: predict whether the current position is a win for the player with the next move (positive) or not (negative). Binary, categorical-only classification. Network: one hidden layer with 64 neurons plus a sigmoid output; 25 epochs, batch size 32.

Vehicle Silhouettes Dataset

(OpenML: vehicle, ID: 54):

846 instances, 18 shape descriptors extracted from 2-D vehicle silhouettes (compactness, circularity, rectangularity, …). Target has 4 classes: bus, opel, saab, van. A multiclass task with numeric features. Network: two hidden layers with 128 and 256 neurons; 4-way soft-max output; 25 epochs, batch size 32.

Vertebra Column Dataset

(OpenML: vertebral-column-2clases, ID: 1524):

310 spine examinations, each described by 6 biomechanical angles (pelvic incidence, pelvic tilt, lumbar lordosis angle, sacral slope, pelvic radius, and spondylolisthesis grade). The target labels a record as normal vs. abnormal, yielding a binary classification task. Network: one hidden layer with 32 neurons followed by a single sigmoid output; trained for 30 epochs with mini-batches of 16.

Wine Recognition Dataset

(OpenML: wine, ID: 187):

178 wines from the same Italian region, each analysed for 13 chemical constituents (alcohol, malic acid, ash, alcalinity, magnesium, …). Goal: assign the sample to one of 3 grape cultivars. Small, all-numeric, 3-class problem. Network: two hidden layers with 64 and 128 neurons; 3-way soft-max output; 30 epochs, batch size 32.

3.2. Classification Results

Table 1, Table 2, Table 3, Table 4, Table 5, Table 6, Table 7, Table 8, Table 9 and Table 10 reveal a strong advantage for fractal activations.3 A fractal function enters the top five of every benchmark and secures the highest mean accuracy on seven of the ten datasets (climate, diabetes, glass, ionosphere, iris, seeds, wine). On the three remaining tasks—vehicle, tic-tac-toe, vertebrae—the best classical activation (ReLU or tanh) leads by at most $1.6$ pp, underscoring the general competitiveness of the new functions.

Dominant Fractal Variants

The modified Weierstrass–Tanh activation is the most reliable performer. It tops four datasets and keeps its standard deviation below $0.04$ in those cases, indicating stable optimisation. Decaying Cosine excels on smaller, low-dimensional problems (glass, iris, seeds) with similar variance. Blancmange is more volatile: under Adagrad on glass the minimum accuracy drops to $0.40$, yet under Adam on wine it attains a perfect $1.00$, suggesting a high payoff for an appropriate optimiser.

Classical Baselines

ReLU remains hard to beat on vehicle silhouettes and in the categorical tic-tac-toe table, while tanh under SGD produces the best vertebrae result. These cases share pronounced class imbalance or geometric structure, where the simplicity of a smooth or piece-wise linear non-linearity may still be advantageous.

Variance and Optimiser Effects

Most large variances trace back to a mismatch between activation irregularity and optimiser dynamics. The modified Weierstrass–Tanh function paired with Adam keeps Std near $0.03$ on seeds and wine, yet with Adagrad on tic-tac-toe its variance triples. Across all tables, 80 % of first-place fractal results use an adaptive optimiser (Adam, RMSProp, AdaDelta), supporting the view that per-parameter learning-rate control is crucial for handling the steeper local slopes of Weierstrass-type curves. SGD, by contrast, performs best when the activation is either smooth (tanh) or only mildly irregular (modified Weierstrass–Tanh), as seen on climate and vertebrae.

Sanity Check: Classical Weierstrass

To rule out the possibility that any “exotic” fractal curve would improve performance by chance, we also evaluated the un-modulated Weierstrass function described in §Section 2.3.5. Unlike its gated or blended relatives, this baseline offers no linear trend and no decay envelope. Across all ten datasets it never enters the top-ten configurations and is usually the worst performer, confirming that fractal structure alone is insufficient. The extra modulation terms—tanh, ReLU, exponential decay—are essential for turning the raw Weierstrass series into a useful activation.

Summary

Fractal activations equal or surpass classical choices on most benchmarks without inflating variance when combined with a suitable adaptive optimiser. Gains are clearest on moderate-size tabular data and signal-rich tasks; strongly imbalanced multiclass problems still favour the economical ReLU. These findings align with the trajectory-length analysis in Section 4: the extra oscillations of fractal curves translate into measurable accuracy improvements provided the optimiser can accommodate their sharper gradients. However, at this point, we acknowledge that our averaged results exhibit variability, i.e., the standard deviation of our averaged results. Still, our findings indicate that fractal activation functions can be a viable choice, as they provide competitive results even within the variability. Furthermore, although there is some variability in the averages, we often achieve the highest maximum accuracies with fractal activations. For some datasets, the variability range even excludes results below the top five, further highlighting their potential.

4. Expressivity of Fractal Activations

To understand why fractal activations help, we measure how much they deform input trajectories inside an untrained, random network. Our protocol follows the length–based diagnostics of Poole et al. [6], Raghu et al. [7]. The experiment has two steps.

1. Generate a probe curve.

Pick two random vectors $x_0,x_1\sim \mathcal{N}(0,I_d)$ and interpolate along the half-circle

$$x\left(t\right)=cos\left(\frac{\pi t}{2}\right)x_0+sin\left(\frac{\pi t}{2}\right)x_1,t\in [0,1].$$

The construction keeps ${\parallel x\left(t\right)\parallel }_2$ essentially constant, so any length change later on stems from the network, not the input.

2. Track length and shape through depth.

For each activation $\varphi$ we initialise an identical fully- connected network ($L=8$ layers, 400 units each, Glorot weights, zero bias). Let $z^{(\ell )}\left(t\right)$ denote the image of $x\left(t\right)$ after ℓ hidden layers. We record two diagnostics:

(a)

Trajectory length${\mathcal{L}}^{(\ell )}={\sum }_{i=1}^{S-1}{∥z^{(\ell )}\left(t_{i+1}\right)-z^{(\ell )}\left(t_i\right)∥}_2,$$S=400$ points. Exponential growth with ℓ indicates that the network is bending and stretching the curve, hence can represent highly non-linear functions.

(b)

PCA (Principal Component Analysis) strip We project ${\left\{z^{(\ell )}\left(t_i\right)\right\}}_i$ onto the first two principal components for layers $\ell =0,1,3,5,8$ to visualise how folds accumulate.

4.1. Observed Behaviour

Table 11 and Fig. Figure 2 and Figure 3 summarise how each activation reshapes a one–dimensional input curve as it propagates through depth $L=8$.

Arc-length growth.

Classical smooth activations behave as expected. Sigmoid shows fast length decay ($8.26\to 1.91$), while tanh grows moderately ($2.8\times {10}^2\to 2.0\times {10}^3$) before saturating. ReLU amplifies length two orders of magnitude faster, reaching $1.1\times {10}^6$ at layer 8. All fractal variants meet or exceed these baselines. The mildest (decaying cosine, modified Weierstrass–Tanh) track the ReLU curve but stop around ${10}^4$. More irregular functions (Weierstrass–Mandelbrot xsinsquared, xpsin, relupsin) explode to ${10}^9$–${10}^{11}$, confirming the strong local stretching predicted by their high-frequency components.

PCA Strips

The projected trajectories in Fig. Figure 2 and Figure 3 complement these numbers. Smooth activations yield gently bent curves; ReLU introduces sharp kinks; fractal ReLU and Weierstrass variants add dense folds that repeat at multiple scales. The pattern is most pronounced for modified Weierstrass–ReLU: the strip shows large global swings overlaid with fine oscillations, reflecting simultaneous expansion and local twisting.

Taken together, the trajectory lengths and PCA strips show that fractal activations not only maintain but extend the expressive capacity of standard networks, provided the optimiser can cope with the steeper local gradients.

Summary

Fractal activations keep the exponential rise in trajectory length found in standard networks, while the PCA plots reveal additional local folds introduced by their oscillatory shape. Networks can therefore access a larger variety of functions without sacrificing depth efficiency.

5. Discussion/Conclusion

Fractal activation functions based on Blancmange and

Weierstrass–Mandelbrot series expand the functional reach of ordinary feed-forward networks. Across ten tabular benchmarks they achieve top-five status in every case and lead seven, with only modest variance when paired with adaptive optimisers. We acknowledge that the averaged results come with some variability, but the results still show that the fractal activation functions are a valid choice that can achieve top performance and/or rank among other well-performing activation functions on average.

Trajectory-length and PCA diagnostics confirm that the new functions preserve the exponential growth in feature complexity associated with depth while adding multi-scale folds absent from ReLU and tanh. Together, these empirical and geometric results show that engineered self-similarity is a viable new design paradigm for neural architectures. Further, our neural network architecture was initially optimized to work well with ReLU activation functions, yet our fractal activation functions still performed best in most instances. This further emphasizes the potential capabilities of fractal activation functions.

Limitations

Our approach is limited to tabular data of rather small sizes adn correspodning shallow neural networks. We showed that fractal activation functions can provide an edge to regular activation fucntions on such dtaa sets wiht these neural network archtiectures. Thus we cnanot infer how these activaiton fucntions would fare on larger data sets. Particualry, ofr e.g. lagnauge proceissing (tranformers) or image processign and classficaiton thereof (vonvolutional neurla networks).

Limitations

Our approach is limited to tabular data of rather small size and the corresponding use of shallow neural networks. We showed that fractal activation functions can provide an edge over regular activation functions on such datasets with these neural network architectures. However, we cannot infer how these activation functions would perform on larger datasets. In particular, this includes applications in language processing (e.g., transformers) or image processing and classification (e.g., convolutional neural networks). While our results suggest increased expressivity for small architectures, and thus the potential to reduce deeper neural networks and make them more efficient for particular use cases, we did not explore deeper networks, leaving the question of how fractal activations scale with depth open. This is a task beyond the scope of the present introduction, but given our results, and particularly the evidence that these activation functions do work for neural networks, it represents the logical next step for our research. We also used fixed fractal hyperparameters (e.g., truncation depth, weights a, b, decay $\lambda$), and did not include a full ablation study on their effects in neural networks, which is a considerable effort beyond the scope of this article.

Future Work

(a)

Depth replacement. Measure whether a shallow network equipped with fractal activations can match or outperform a deeper classical model on the same task. Such a study would clarify how much “effective depth’’ the fractal folding provides.

(b)

Fractional optimisers. Explore gradient methods that use fractional derivatives of the loss surface [29,30,31,32,33]. Their built-in memory and/or regularization could counteract the potential large local slopes that arise from fractal activations and further stabilize training.

(c)

Parameter search. Systematically tune the geometric weights $(a,b)$ of the Weierstrass sums and the envelope decay $\lambda$ to balance expressivity and gradient magnitude. A Bayesian search over this small hyper-space may yield task-specific defaults.

(d)

Complexity/Multidimensionality of a Dataset Another interesting perspective would be to analyze datasets for their complexity, e.g., using support vector machines and counting the kernel dimensions, and then relating this sort of complexity/data dimensionality to the usage of fractal activation functions for normalized neural network architectures. I.e., finding out whether the induced fractal folding can compensate for high-dimensional data complexity and thus improve accuracy.

(e)

Convolutional and/or Transformer Architectures

Addressing these points will clarify when fractal activations are most useful and how they can be combined with optimisers that exploit the same long-range dependencies that define the functions themselves.