Multifractal Properties of BK Channel Currents in Human Glioblastoma Cells

Potassium channels play an important physiological role in glioma cells. In particular, voltage- and Ca2+-activated large-conductance BK channels (gBK in gliomas) are involved in the intensive growth and extensive migrating behavior of the mentioned tumor cells; thus, they may be considered as a drug target for the therapeutic treatment of glioblastoma. To enable appropriate drug design, molecular mechanisms of gBK channel activation by diverse stimuli should be unraveled as well as the way that the specific conformational states of the channel relate to its functional properties (conducting/nonconducting). There is an open debate about the actual mechanism of BK channel gating, including the question of how the channel proteins undergo a range of conformational transitions when they flicker between nonconducting (functionally closed) and conducting (open) states. The details of channel conformational diffusion ought to have its representation in the properties of the experimental signal that describes the ion-channel activity. Nonlinear methods of analysis of experimental nonstationary series can be useful for observing the changes in the number of channel substates available from geometrical and energetic points of view at given external conditions. In this work, we analyze whether the multifractal properties of the activity of glioblastoma BK channels depend on membrane potential, and which states, conducting or nonconducting, affect the total signal to a larger extent. With this aim, we carried out patch-clamp experiments at different levels of membrane hyper- and depolarization. The obtained time series of single channel currents were analyzed using the multifractal detrended fluctuation analysis (MFDFA) method in a standard form and incorporating focus-based multifractal (FMF) formalism. Thus, we show the applicability of a modified MFDFA technique in the analysis of an experimental patch-clamp time series. The obtained results suggest that membrane potential strongly affects the conformational space of the gBK channel proteins and the considered process has nonlinear multifractal characteristics. These properties are the inherent features of the analyzed signals due to the fact that the main tendencies vanish after shuffling the data.


■ INTRODUCTION
The available chemotherapy and radiology treatments turn out to be ineffective in the case of gliomas, which are brain tumors arising from glial cells. 1,2 Gliomas account for the majority of malignant brain tumors in adults 3 and are graded from I to IV, with higher grades being more differentiated and malignant. 2 Grade IV gliomas are called glioblastoma multiforme (GBM) and exhibit the highest proliferative potential and almost complete resistance to currently available therapies. 4,5 The great majority of patients with GBM do not survive beyond 2 years even when a combination of novel and conventional therapies, i.e., surgery, chemotherapy, and radiotherapy, was introduced. 1,4−6 Focal surgical resection is ineffective and adequate radiotherapy is impossible in glioblastoma due to the fact that GBM is characterized by extensive invasion, migration, and angiogenesis. 7 To enable the development of a more effective therapeutic treatment for glioblastoma, one should better understand all biological processes taking place at the molecular and cellular levels in this tumor. Several ion channels have been implicated in glioblastoma proliferation, migration, and invasion. 7 In this work, we focus on the activity of voltage-and Ca 2+activated large-conductance K + channels (BK) obtained from human glioblastoma cells (gBK channels). BK channels are overexpressed in malignant gliomas (in comparison to nonmalignant cortical tissues), and their expression level correlates positively with the malignancy grade of the tumor. 8−10 These channels are expressed in specific isoforms in glioma cells that have slightly different characteristics from other BK channel exons (e.g., they are more sensitive to cytosolic concentration of calcium ions 11 ); thus, they are called gBK channels. The gBK channels play an important physiological role in glioma cells: they program and drive cell growth and extensive migration (also in glioblastoma stemlike cells), 7,9,10,12,13 so they detrimentally facilitate the invasiveness of glioblastoma that renders them incurable so far.
One can indicate several processes taking place at the molecular level, in which gBK channels can contribute to the shape and volume changes of glioma cells during their invasive migration in a crowded environment (by triggering and motorizing that process). 12 First, the distribution of ions affects the water flow across the cell membrane, and consequently, the cell volume, due to the effective osmotic pressure. 14 Taking into account the overexpression of gBK channels in gliomas, these channels can detrimentally affect osmosis and regulation of cell volume. Second, gBK channels are anticipated to provide the electrochemical driving force for the ion movement needed for the release of cytoplasmic water and cell shrinkage, which in turn facilitates the extensive migrating behavior of glioblastoma cells. This was indicated by the results provided in the works of Wondergem et al. 15,16 Moreover, due to the mechanosensitivity of gBK channels, they may be involved in the mechanotransduction during cell shape and volume changes. 17 Also deformation or reorganization of the cytoskeleton during an alteration in cell volume or shape may affect the functioning of gBK channels, 17 which, as a consequence, should affect the effectiveness of cell migration.
Here, we go beyond the usual description of the changes of gating kinetics resulting from the electrical stimulation of cell membrane and investigate the effects of membrane depolarization on a channel's conformational dynamics. Mechanisms of K + channel activation by diverse stimuli like voltage, Ca 2+ , Mg 2+ , H + , HEME, temperature, mechanical strain, etc. still evoke open discussions among researchers. The authors who model activation dynamics and gating present many possible scenarios, which differ by the number of available channel states within both conducting (functionally open) and nonconducting (functionally closed) states' manifolds; also the mechanics of the state transitions are introduced in the models in a different way. 18−24 A great breakthrough has been made by resolving the molecular structure of the Ca 2+ -bound and Ca 2+ -free BK channels. 25−27 It allowed not only for the inference that voltage and Ca 2+ sensors are coupled and they can cooperatively influence the channel pore gate domain, but also these studies failed to identify a physical gate that could mechanically block the ion flow in the nonconducting state. 28 From this point of view, the pore-forming helices do not form a tight bundle as in some Shaker-like channels to dam K + transport through the pore. In ref 25, the authors observed four different structures for the BK channel when the channel was neither voltage-activated nor Ca 2+ -activated. In those terms, the channel would be expected to be functionally closed. Nevertheless, the inner poreforming helices are not so close to each other to prevent potassium ions to access the selectivity filter. Quite similar observations are made in the studies on the activation of ligandgated K + channel from the Slo family. 29 Namely, the Na +dependent K + channel Slo2. 2 25,26,29,30 one can expect that at least some of the aforementioned inferences may also refer to the BK channel's conformational dynamics.
BK channels pass through multiple kinetic states over time during gating. It can be assumed that Hite et al. 25,29 have identified some of the structures corresponding to different kinetic states, providing insight into which conformational states' functionally closed and open states might be adopted. 28 Still, some questions arise, among others: • Where is the actual channel gate that allows for conducting/nonconducting fluctuations at fixed conditions? • What kind of and how many stable structures of the BK channel protein exist at intermediate Ca 2+ and/or voltage levels? • Are there any differences in system dynamics in functionally open and closed states, and which ones influence the system as a whole to a higher extent? A possible answer for the first question is provided by the hydrophobic gating mechanism postulated in ref 31, where the authors conclude that the BK channel does not need a physical gate. Spontaneous flickering between conducting and nonconducting states at fixed conditions can be realized as wetting and dewetting of the channel pore resulting from the fact that the pore can constantly undergo changes in shape and surface hydrophobicity (conformational diffusion). The answers for the above questions need broader studies and discussion. However, some clues may be provided here by means of nonlinear analysis of experimental time series describing single channel activity.
As shown in our previous research, 17 application of nonlinear methods in the analysis of an experimental nonstationary series describing ion-channel activity can be useful for estimating the changes in the number of channel states available from a geometric and energetic point of view at given external conditions. We claim that channel conformational diffusion has its representation in the structure of the signal. The most structurally distant conformations are suspected to affect the complexity of the experimental data at least. In this paper, we would like to continue and broaden the discussion about the conformational diffusion of the gBK channel in electrically stimulated membranes based on the results of the multifractal detrended fluctuation analysis (MFDFA). 32 For over 30 years, ion-channel recordings have been analyzed paying special attention to their nonlinear, fractal properties. 33−38 Here, we not only describe the multifractal nature of the experimental time series of BK channel activity but also present a biological interpretation of the obtained characteristics, which gives an additional insight into the conformational dynamics of the channel protein.

■ METHODS
Detrended Fluctuation Analysis (DFA). The DFA method was first proposed by Peng in 1994 for investigating the correlation in DNA structure. 39 The last years have seen a renewed importance in the application of this method to biological data and also as a technique capable of distinguishing between healthy subjects and heart failure patients. 40 This technique relies on the assumption that the signal is influenced by both short-term and long-term features. For a proper interpretation of the effects hidden behind the internal dynamics, the signal is supposed to be analyzed at multiple The Journal of Physical Chemistry B pubs.acs.org/JPCB Article time scales. 41 A brief description of the original DFA algorithm is presented below. The procedure starts with the calculation of the profile y i as the cumulative sum of the data x i with the subtracted mean ⟨x ⟩ Next, the cumulative signal y i is split into n s equal nonoverlapping segments of size s, for which we use powers of 2, s = 2 r , r = 4,..., 11. For all segments of size s, v = 1,..., n s , and the local trend y v,i m is calculated. In the standard DFA method, the local trend is calculated by means of the least-squares fit of order m. In this work, a second-order polynomial m = 2 is used. The variance F 2 (s, v) as a function of the segment length s is calculated for each segment v separately.
As the last step, the Hurst exponent H can be determined as the slope of the regression line of the double-logarithmic dependence log F(s) ∝ H log s of the square root of the average variance A similar power law S(q, s) ∼ s H(q) can be utilized to calculate the q-order generalized Hurst exponent H(q). The latter is required to compute the singularity spectrum for a time series. In the first step, the mass exponent τ is determined via the relation τ(q) = qH(q) − 1. Next, the qth-order Holder exponent h(q), a quantity that characterizes the singularities, is estimated as a derivative of the mass exponent, h q ( ) Finally, the qthorder multifractal singularity spectrum D(h) (mf-spectrum) can be constructed as a Legendre transform of the mass exponent It results in the usual concave-shaped distribution of the singularity strengths. 43 Focus-Based Multifractal Formalism. Multifractal properties of a time series are reflected in the generalized Hurst exponent H(q). The standard method for the estimation of those properties often results in a nonmonotonic characteristic of this function, which causes the degeneracy of a singularity spectrum. 44,45 It is typically a consequence of the finite length of the series. The scaling function S constructed from different moments q of the measure depends on the scale of observationsee eq 3. However, if one considers the full length of the time series, the dependency on the moment q disappears. It means that regardless of the moment, there exists only one value for the fluctuation function S(q, L) = S(L) = S L . It defines a single pair of values (log S L , log s) to which the function log S(q, s) for any moment q supposes to converge with a growing scale s. We will call this pair of values the focus point. The existence of this point is the central idea of the focus-based multifractal (FMF) formalism. 46 The traditional approach to find a singularity spectrum is based on finding a linear representation of the calculated log S versus log s characteristics for each moment q. An algorithm, however, does not guarantee that the such-determined q− dependent linear functions will cross with each other at the focus point or even will not cross at all (see the blue solid lines in Figure 1 for details). Instead of calculating separate fits for selected moments q at a time, one can try to fit the whole family of functions S(q, s) keeping the requirement of sustaining one focus point (see Figure 1). The usual way is to construct a cost (loss) function that would capture the properties of both the focus point and the whole family of the qth-order scaling functions. This function can be deduced from the power law If we mark an iteratively updated scaling function with a hat, then the total mean squared error cost function MSE to be minimized reads 46 where n s is the total number of segments of size s and n q stands for the number of scaling factors q. Next, the scaling functions are found for the data. In the subsequent step, the iterative optimization algorithm, e.g., the first-order gradient descent, is utilized to find the minimum of the cost function.
The newly obtained family of the linear representation of the scaling functions can be later used to calculate a new generalized moment-wise Hurst exponent function Ĥ(q), which in turn directly yields the singular spectrum D(h). From this point onward, the procedure for the calculation of the singular spectrum is the same as that for the classic MFDFA algorithm. 42 The analysis by means of the standard MFDFA and focus-based method was performed using our own authorial Python software, where we implemented the methodology developed by Mukli et al. 46 The parameters of the singularity spectra (viz., spectral width Δ, half-width Δ 1/2 , maximum of a spectrum H max , spectral symmetry) obtained for the experimental time series of single channel currents describe the complexity and multifractal properties of these data. Such characteristics can be, in turn, interpreted in terms of the process underlying the observed current fluctuations, which is the conformational dynamics of a channel protein. An appropriate discussion will be provided in the Results and Discussion section.   The Journal of Physical Chemistry B pubs.acs.org/JPCB Article 2, the probability of conducting state (O) increases with membrane depolarization, which is typical for voltage-dependent channels. Single-channel current amplitude increases as the difference in membrane potential increases; thus, local maxima in Figure 3 become well separated at both deep depolarization and hyperpolarization. Event Detection. The MFDFA is carried out on a raw experimental data seriestime series of single-channel currents and also on preprocessed resultsseries of ionic currents recorded during the conducting (open) state of a channel and some current fluctuations recorded during the nonconducting state of a channel. The threshold current value used to identify transitions between the subsequent states is evaluated as given in ref 36. The analysis of single-state data sets was performed to indicate the complexity of the signal within currents corresponding to a single manifold of states and any differences between the properties of the conducting and nonconducting states, as well as to determine which kind of signal corresponding to conducting or nonconducting states determines the characteristics of the total experimental time series of single-channel currents.

■ RESULTS AND DISCUSSION
The necessity of applying a modified MFDFA technique in the form of an FMF analysis is presented in Figure 4. The obtained spectra were calculated by both the methods: standard MFDFA (blue dots) and the focus-based modification (orange crosses) (see Figure 4). In the standard method, the degeneration of the spectrum (zigzag type) is found. In other words, a large number of cases are characterized by the inversed singularity spectrum. The nonmonotonic relationship of the generalized Hurst exponent as a function of q precedes obtaining this type of spectrum shape. This situation causes difficulties in the interpretation of multifractality of the data and also the proper estimation of spectrum parameters such as maximum of spectrum or the spectrum width.
Comparison of mf-Spectra at Different Membrane Potentials. A comparison of mf-spectra calculated by the FMF method is presented in Figure 5. For the investigation of the long time series (500 000 data points), the range of scales s ⊂ [2 4 ,2 13 ] was selected. The comprehensive analysis of the nature of the BK channel's multifractality requires the calculation of the spectra for the following cases: (i) raw data and (ii) data after the shuffling operation. Figure 5 summarizes the spectra obtained for all subsequent channel states (total signal composed of the experimentally recorded ionic currents), series of potassium currents corresponding to the conducting state of a channel, and series of some "leak" currents that correspond to the nonconducting states of a channel at different stages of depolarization and hyperpolarization. During the analysis of the total signal, one can note a marked tendency between the different values of membrane potentials in both cases, at hyperpolarization and depolarization. The average spectra of the data obtained when the value of membrane potential was closest to zero in each group (20 and −20 mV) are clearly shifted to the smaller values of h(q). In other words, the spectral maximum is most extended to left at this specific condition, and then along with increasing applied potential at membrane depolarization and decreasing applied potential at membrane depolarization successively moves toward larger values of h(q) (Figure 5a,b). Considering the results of the multifractal analysis dedicated to the currents recorded during the open (conducting) and closed (nonconducting) states of the channel, one can observe that the results obtained during the channel's closures are completely consistent with those corresponding to the total signals.
The results characterizing open states are the opposite of the remaining ones. First, the variability of their spectral width is substantially smaller. Second, the general trend shows an increase in spectral width when the membrane potential decreases (only for −40 mV, there is a local minimum). Such results suggest that the total signal's characteristics are mainly determined by the recordings obtained during the nonconducting states of a channel. The recognized differences in spectral distributions allow us to infer that the dynamics of conformational changes within the conducting and nonconducting states' manifolds differ significantly. Roughly speaking, single-channel currents recorded when the channel pore exhibits possibly high conductance retain a self-similar structure over a range of scales regardless of the membrane potential ( Figure 5b) (which is also notable by the trendreinforcing behavior, as measured by Hurst exponents). Whereas current fluctuations recorded when the channel is supposed to not conduct potassium ions as well as the total signal significantly lose their multifractal self-similarity (and the recorded time series become uncorrelated or even anticorrelated, as shown by values of the Hurst exponents) (Figure 5a,c).
It is also visible in the case of functionally open states that multifractal spectra have a right truncation. A long left tail suggests that the time series of channel currents have a multifractal structure that is insensitive to the local fluctuations with small magnitudes.
After the shuffling operation (Figure 6), the spectral distribution is quite different, but the most important aspect here is that the spectra are about twice narrower. The effect of the reduced multifractality after mixing of the data has a large consequence in the proper interpretation of the source of the fractal nature of the examined signals.
There exist two general sources of multifractality which have influence on the shape of the mf-spectrum: (i) the broad probability density function (pdf), which lies behind the data, and (ii) different behaviors of the (auto)correlation function for large and small fluctuations. Furthermore, both situations are possible simultaneously. In case (i), shuffling will not change the mf-spectrum; for (ii), it will destroy the effect completely as the Figure 4. Comparison of mf-spectra for standard multifractal formalism (blue dots) and modified focus-based method (orange crosses) calculated for the representative measurement obtained at membrane hyperpolarization, U = −20 mV. The black and red dots mark, respectively, the maximum of the spectrum and the generalized Hurst exponent for the nondegenerated case.
The Journal of Physical Chemistry B pubs.acs.org/JPCB Article shuffling will erase the possible correlations. When cases (i) and (ii) are present simultaneously, the spectrum will differ from the original one, and the weaker multifractality can be identified. In our case, we can observe exactly the last mixed situation, and thus we suspect that the multifractality of the data is caused by both the correlation and the broad pdf.
To interpret these results in terms of the considered biological system, one has to note the following facts: • BK channels are voltage-activated, which means they exhibit more often the conducting state than the functionally closed one at membrane depolarization and tend to retain a nonconducting state at hyperpolarization of the cell membrane. Nevertheless, even at a negative potential, there is a nonzero probability that the channel rapidly opens for a relatively short time (as shown in Figure 2), • the actual single-channel conductance is determined by two factors. The first, and the most detrimental effect, is exerted by the membrane potential. The absolute value of a single-channel current increases as the difference of the electric potential on both sides of the channel membrane increases (both toward highly positive values and toward negative ones), as shown in Figure 2a. Higher amplitudes of channel currents at deep-membrane depolarization or hyperpolarization result in broader probability density functions (Figure 2b). Second, the conductance of the channel pore in the open state varies with voltage as a result of the structural changes that a given channel undergoes during voltage activation. 25,26,47,48 But these slight changes in geometry can significantly influence the kinetics of switching between conducting and nonconducting states.
Taking into consideration the facts mentioned above, multifractal properties of the analyzed time series of channel currents (total signal) are an inherent feature of the system connected to the dynamics of switching between the channel states, which change with membrane potential but does not depend strictly on the value of current amplitude. The changes in channel currents with voltage pertain mainly to the conducting states of a channel since the nonconducting states form a baseline during the experimental recording and they underlie smaller fluctuations in the examined range of the membrane. Single-channel currents The Journal of Physical Chemistry B pubs.acs.org/JPCB Article recorded at the conducting states of the channel have multifractal characteristics over a range of scales regardless of the membrane potential ( Figure 5b); in contrast, the multifractality within closed states and the total signal vary significantly with voltage. As the picture of multifractal properties of the total signal fits to the one obtained for functionally closed states, one can infer that the channel dynamics is mainly influenced by the dynamics of conformational transitions within the nonconducting states. This inference is compatible with some popular kinetic models of the ion-channel activity, 22,49,50 where more kinetic substates correspond to functional closures than openings of a channel. It is possible that the scheme of switching between functionally closed conformations becomes more complex with the increase of the difference in electric potential on both sides of the membrane, which leads to an eventual widening of the mfspectrum at these conditions. The confirmed existence of the second source of multifractality of the investigated data, namely, correlations for large and small fluctuations, is also worth noting. The analyzed time series are long-term correlated at all experimental conditions for channel currents at the conducting states of a channel and at high differences of electric potential on both sides of the membranefor nonconducting states, and total signal, as shown by the values of Hurst exponent (Figures 5 and 6).
A strong analogy exists between the multifractal and thermodynamical characteristics. In particular, multifractal spectra can be related to entropy. 51−53 In total signal as well as in the series of channel currents in functionally closed states, one can observe shifting of the maximum of multifractal spectra toward higher values when the difference in electric potential on both sides of the membrane increases (both at depolarization and hyperpolarization). It indicates greater complexity of the signal at highly positive and negative potentials comparing with the data obtained at membrane potentials close to zero, which may suggest an increase in the number of attainable channel substates (mainly within the nonconducting manifold) with absolute value of voltage. The symmetry of the changes in signal multifractality (and, consequently, entropy) occurring both at membrane depolarization and hyperpolarization may be counterintuitive in the case of a voltage-activated channel. The findings from refs 22, 49, 50, and 54 suggest nonsymmetric nets of conducting and nonconducting states, but there is no information about probabilities of switching between different

■ CONCLUSIONS
In this work, a novel approach of multifractal signal analysis is presented. To the authors' best knowledge, very few publications can be found that discuss the issue of the multifractal character of a time series by implementing an FMF methodology. An implementation of this unique technique, which is capable of handling empirical signals with a varying degree of heterogeneity, brings a lot of valuable information to the investigation of an ion channel's activity. This work concludes that the multifractality can be regarded as an inherent feature of the single-channel currents obtained by patch-clamp measurements. The observed multifractal spectra suggest that the characteristics of system dynamics are substantially different in functionally open and closed states, and the total signal recorded during experiments is influenced to a higher extent by the nonconducting states than the conducting ones. It is quite interesting that the symmetric increase of spectrum width and shifting of maximum of mf-spectrum of both the total signal and channel currents recorded during the functionally closed states toward higher values as the difference in electric potential on both sides of the membrane patch increases. It suggests a higher complexity and entropy of the signal recorded at both strong membrane depolarization and hyperpolarization comparing with the ones obtained at moderate membrane potentials. According to Boltzmann's definition of entropy, these results ought to indicate an increase of the attainable substates (stable conformations) mainly in the nonconducting domain with the absolute value of applied voltage. Regardless of the applied voltage, the time series of channel currents recorded at the channel's conducting conformations are nonrandom but caused by the orderly process exhibiting long-range correlation features. To sum up, channel dynamics are qualitatively and quantitatively different in the case of conducting and nonconducting states of a channel. Assuming that the most distant conformational states from energetic point of view have their representation in the recorded signal (single-channel current), it can be noted that the states obstructing ionic flow through a channel pore are more complex and influence the multifractality of the total signal to a higher extent than the ones allowing for K + transport. One should remember that our analysis does not discern between mechanically closed and nonconducting open statesso both groups, physically blocked conformations and sufficiently narrow ones (implying hydrophobic gating), predominate in shaping the channel's activity patterns. An interesting task for future investigation can be to carry out a comparative MFDFA analysis of a patch-clamp time series of a single-channel current on a wild-type and genetically modified BK channel that cannot exhibit relatively narrow conformations of the pore enabling for hydrophobic gating. Such analysis could be used to discriminate between the impacts of both the aforementioned groups of nonconducting states on the total signal. Moreover, the presented multifractal analysis can be a tool of supplemental analysis procedures, e.g., in cases when one should determine to what extent different regulatory β subunits can modulate the complexity of channel behavior. The results of such analysis could help answer the question whether they modulate the relative stabilities of preexisting conformations 27 or create new ones.
In cases of glioblastoma, current medical approaches turn out to be almost powerless. Among the challenges to curing primary brain tumors, one can list the development of a precision medicine approach to treating brain tumors. In that aspect, novel approaches should be introduced, which could be based on artificial intelligence (AI) (e.g., deep learning, neural networks). The AI methods can be used for diagnosing, managing, and designing drugs against gliomas. In the literature, there already exist some reports like ref 55, where the authors present novel AI approaches to predict the grading and genomics from imaging, automate the diagnosis from histopathology, and provide insight into prognosis. Taking into account the therapeutic potential of gBK channel modulators in the treatment of glioblastoma, one could propose some AI methods to determine a group of active substances that could act as a drug against gliomas. Machine learning might be developed with the aim to determine patterns within the experimental data describing ion-channel activity, where some of the classification algorithms could be based on the results of an MFDFA analysis. (Our preliminary analyses suggest that multifractal analysis better discriminates singlechannel current from different exons of the BK channel than do the kinetic characteristics.) Complexity and multifractality of a signal describing different BK channel exons bound or unbound to ligand molecules (specific modulators) could be one of the factors used in optimizing the structures of potential BK channel modulators used as a drug against glioblastoma.