Differences in Gating Dynamics of BK Channels in Cellular and Mitochondrial Membranes from Human Glioblastoma Cells Unraveled by Short- and Long-Range Correlations Analysis

The large-conductance voltage- and Ca2+-activated K+ channels (BK) are encoded in humans by the Kcnma1 gene. Nevertheless, BK channel isoforms in different locations can exhibit functional heterogeneity mainly due to the alternative splicing during the Kcnma1 gene transcription. Here, we would like to examine the existence of dynamic diversity of BK channels from the inner mitochondrial and cellular membrane from human glioblastoma (U-87 MG). Not only the standard characteristics of the spontaneous switching between the functional states of the channel is discussed, but we put a special emphasis on the presence and strength of correlations within the signal describing the single-channel activity. The considered short- and long-range memory effects are here analyzed as they can be interpreted in terms of the complexity of the switching mechanism between stable conformational states of the channel. We calculate the dependencies of mean dwell-times of (conducting/non-conducting) states on the duration of the previous state, Hurst exponents by the rescaled range R/S method and detrended fluctuation analysis (DFA), and use the multifractal extension of the DFA (MFDFA) for the series describing single-channel activity. The obtained results unraveled statistically significant diversity in gating machinery between the mitochondrial and cellular BK channels.

It has been shown that the mitochondrial variants of the BK channels (mitoBK) are expressed when the Kcnma1 gene undergoes splicing to the DEC isoform [19]. The DEC splice variant is recognized by the presence of an insertion of a 50-aa C-terminal sequence named after the three last amino acids, which prevents the expression of the BK channel at the plasma membrane. In so far as the localization of the BK-DEC splice variant of BK-type channels is confirmed in mitochondria, it is still not known whether this isoform imposes functional changes in channel's activity compared to the other BK isoforms present in the cell membrane. Additionally, differences of biophysical properties between plasma membrane-and mitochondrial channels can result from various factors, including the presence of membrane-specific proteins or differences in the lipid composition of the inner mitochondrial membrane in comparison with the plasma membrane. The inner mitochondrial membrane is characterized by a high degree of unsaturation, the presence of cardiolipin, and the absence of cholesterol, whereas the plasma membrane is highly saturated, contains relatively large amounts of sphingomyelin and cholesterol [20]. These discrepancies between the mitochondrial and cellular membranes can also escalate the differences in gating between the BK channel isoforms located in the plasma membrane and their mitochondrial analogs [21][22][23][24][25]. Differences in the single-channel conductance can be observed, as in the case of the ROMK-like channels from the plasma membrane and their analogs in the inner mitochondrial membrane [26].
In this work, we analyze the differences in the single-channel gating dynamics between the BK channels from the inner mitochondrial and cellular membrane from human glioblastoma (astrocytoma, U-87 MG). The former reports [27,28], where the results of Western blot analysis, immuno-gold electron microscopy, high-resolution immunofluorescence assays, and/or polymerase chain reaction were summarized, determine the presence of BK channels in U-87 MG cell line and give strong evidence that only BK-β4 assemblies were present in the analyzed patches of the cellular and mitochondrial membrane. Moreover, the presence of BK channels was also confirmed by the use of paxilline and iberiotoxin which exerted a channel inhibiting effect in appropriate patch-clamp experiments. In turn, the application of BK channel opener NS1619 caused the increase of open state probability up to 0.95. Regulation of these channels by calcium were also confirmed [27].
Due to the fact that the BK channels are expressed in specific isoforms in glioma cells which differ from the other BK channel isoforms by the enhanced sensitivity to cytosolic Ca 2+ concentration, they are called gBK channels [29]. The gBK channels in the plasma membrane play an important role in the physiology of glioma cells. They are involved in cell growth and extensive migration (also in glioblastoma stemlike cells), so they facilitate the invasiveness of glioblastoma [28,[30][31][32][33]. Considering the putative physiological meaning of the mitoBK channels in gliomas, it was reported that they can be involved in the regulation of the respiratory chain [27]. Moreover, the involvement of mitoBK channels in the phenomenon of cytoprotection has been demonstrated in [34][35][36]. In consequence, activity of both groups of BK channels in gliomas is one of the component processes which renders this kind of cancer incurable so far.
From this perspective, the gBK channels seem to be a promising candidate as a drug target in novel therapies against glioblastoma, which is particularly important since even a combination of contemporary therapies (i.e., surgery, chemotherapy, radiotherapy, or even immunotherapy) turns out to be not sufficiently effective in glioblastoma treatment [37][38][39][40]. Development of highly specific modulators acting most effectively on the gBK channels within cell membranes could prevent shape and volume changes of glioma cells which would hamper triggering and motorizing their invasive migration in a crowded environment. Whereas targeting the mitochondrial BK channel variants may help to induce cell death by e.g., hampering its respiration. This kind of highly specific channel-targeted therapy is a challenge for molecular pharmacology. Nevertheless, to reach this milestone goal and avoid the side-effects [41,42], complex investigations should be carried out in a step-by-step manner, where subsequently the details of the whole molecular machinery of gating should be unraveled, and then possible active sites, as well as the highly specific ligands, should be indicated. We are convinced that the first task to realize that aim is to answer the question of whether there exist any differences in gating dynamics of BK channels from mitochondrial and cellular membranes, which we study in this work. The possible dynamical diversity can help to adapt the mitochondrial and cellular channels to the local conditions and ensure the maximal effectiveness of fulfilling their physiological role. In turn, the description of the changes between the patterns of conducting/non-conducting channel fluctuations can bring useful information which can be utilized in resolving the mechanism of switching between channel substates (relatively stable conformations).
The experimental data, which are analyzed here, are obtained by the patch-clamp method in a single channel mode at external conditions adjusted to impose possibly similar open probability levels of the channel for both cases-mitochondrial and cellular patches. Due to the significant differences between the Ca 2+ sensitivity of the plasma membrane and mitochondrial BK channels indicated by us during a preliminary experimental set and confirmed by the studies from the literature [18,29,[43][44][45], we applied different concentrations of Ca 2+ in bath/pipette solutions to cover the open state probability variability range from below 10% to over 90% at the applied membrane potentials from −60 mV to 60 mV with a 20 mV step for both BK channel variants. To obtain a detailed picture of gating dynamics, first, we compare the recorded traces by the construction of voltage-activation curves. Then, we perform the correlation analysis of the time series of single-channel currents (all currents obtained experimentally, and the series of the currents corresponding to functionally open and closed channel separately) and the dwell-times of subsequent channel states (in three variants, as previously, i.e., all states, only conducting and only non-conducting ones). Despite this kind of analysis remains still quite unconventional in the field, we decided to apply it to the experimental data describing single-channel activity because of the large interpretative potential of the calculated characteristics. Following measures of the short-and long-range correlations are evaluated: • the dependencies of mean dwell-times of (conducting/non-conducting) states on the duration of the previous state (here called "conditional mean dwell-times"), • Hurst exponent by the R/S method [46], • generalized Hurst exponent (scaling exponent) by the detrended fluctuation analysis (DFA) [47], • spectral parameters obtained by MFDFA analysis [48].
The obtained results can be interpreted in terms of the complexity of the switching mechanism between stable conformational states of the channel of different pore-geometry and lifespan. Thus, our findings can indicate whether the splicing sequence corresponding to a particular location and/or the location-specific composition of the membrane can define the conformational machinery of the channel gating. The Hurst R/S analysis describes a long-range memory effect in analyzed series [46]. This method was applied to the analysis of the single-channel activity in [49][50][51][52][53][54] where the authors analyzed series of dwell-times of subsequent channel states corresponding to experimental results, or the results of simulations of channel gating models which should reproduce this feature of the real system. It turns out that the analyzed series of channel states' durations are long-range correlated. The short-lasting states tend to be followed by short ones, and long-lasting states -by long ones. This tendency is conserved over the range of time scales. The physiological meaning of this feature, as well as its origin, remain still not understood. A possible explanation of the long-range memory effect in channel dynamics was proposed in [54], where the authors postulate that the single-channel activity is highly biased by another independent process which occurs at a larger time scale in relation to the conformational diffusion of a channel protein, but still, it can deeply affect the diffusive space of the channel gate. A reasonable candidate for such a component process are fluctuations in membrane thickness near the channel's location. Squeezing and relaxation of the membrane ought to detrimentally influence the mass density of the channel protein, and, thus, change the number of energetically equivalent substates within a given macrostate (conducting/non-conducting state).
The investigation of correlation properties of ion channel fluctuations by means of the Detrended Fluctuation Analysis (DFA) was previously described in [55] where the authors analyzed the series of dwell-times of subsequent channel states and in [56] where time series of single-channel currents were analyzed (the raw signal and, separately, series of those current values which correspond only to an open channel state, or only to a closed channel state). Both reports indicate persistent long-range correlations in ion channel gating and stronger correlation effect within the closed states than the conducting substrates of the channel [55,56]. The single-channel currents from the voltage-dependent K + channel in rat dorsal root ganglion neurons revealed the existence of memory effect by means of the DFA analysis [57].
The Multifractal version of the Detrended Fluctuation Analysis (MFDFA), where the basic DFA technique is extended over the range of statistical moments of the calculated variance in terms of the scaling function [48], was also used for the study of ion channel's activity in [58,59]. The reported results addressing either dwell-time series of channel states [58] or single-channel currents [59] indicated that the channel gating is an orderly process exhibiting long-range correlation features of a multifractal form. The membrane potential strongly affects the obtained spectral parameters. The spectrum width corresponding to the experimental signal increases with the rising difference in electric potential on both sides of the membrane patch [59]. Such a result suggests that a higher complexity and entropy of the signal is reached both at strong membrane depolarization and hyperpolarization when compared to the recordings obtained at moderate membrane potentials. This, in turn, may be interpreted in terms of changes in attainable substates (stable conformations) at given experimental conditions, according to Boltzmann's definition of entropy. Thus, one expects an increase in attainable channel's conformations with the rise of the absolute value of membrane potential [59]. The system dynamics are substantially different in functionally open and closed states and the total signal recorded during experiments is affected to a larger extent by the nonconducting states than the conducting ones [59] (which is in agreement with the inferences from [55,56]).
What is worth noticing the aforementioned results and inferences pertaining to the analysis of the BK channels activity in the cell membrane. The Hurst, DFA, and MFDFA analyses were not conducted for the mitochondrial BK channels hitherto. Thus, the current investigations should supplement the literature in this aspect.
To gain possibly broad information about the channel dynamics from the correlation analyses performed in this study, we apply the Hurst R/S and DFA methods in investigations of both: series of channel currents and the corresponding dwell-times of subsequent channel states. These two approaches allow us to gain characteristics of channel gating dynamics which should be interpreted in a slightly different way. The current fluctuations should resemble the subtle changes of protein conformational substrates (at least between the conformations which are sufficiently remote from each other from the energetic and spatial point of view). Thus the correlation analysis of the single-channel currents should describe the memory effect in the process of switching between conformational states of different pore conduction (i.e., between those conformations that imply different pore geometry, and in consequence, different ability to transport K + ions). For persistent correlation, the high values of the current tend to be followed by high values in the consecutive time step and low values of the current tend to be followed by low values. In turn, during the R/S and DFA analyses of the dwell-time series the main focus of attention is on the character of the transitions between conducting and non-conducting states of different lifespan. The obtained correlation measures characterize the memory effect in switching within the spectrum of channel protein substrates occurring on a wide time scale range. In that case, the possible presence of long-range memory suggests that the short-lived substrates tend to be followed by the short-lived ones, and the long-lived substrates tend to be followed by the long-lived ones.
The former reports in the literature suggest that the memory effects may be observed both within the whole experimental series describing the channel's fluctuations between the subsequent conducting and non-conducting states, and also independently within the functionally open and closed states separately [55][56][57][58][59]. Thus, in this work we analyze the series of channel currents and also the series of dwell-times of subsequent channel states in three variants: the total activity of the channel (which is a superposition of two alternate signals describing closed and open state dynamics), as well as the series where the original traces were divided into two signals, containing only closed or only open state fluctuations, respectively. In this way, we ought to get the information on whether the movements of functional domains within the channel protein that cause gating are correlated regardless of whether the channel is open or closed.
In the case of the MFDFA methodology, we only consider the raw experimental patch-clamp recordings due to the requirement for the high length of the analyzed series. Nevertheless, this kind of analysis has a valuable interpretive potential, since the changes in spectral width refer to the changes in attainable substrates of a channel (stable conformations) at given experimental conditions [59].

Cell Culture
The cell culture was carried out on a human glioblastoma cell line (astrocytoma U-87 MG) as given in [27,59]. Briefly, the cells were cultured in DMEM supplemented with 10% fetal bovine serum at 37 • C in a humidified atmosphere with 5% CO 2 . The other substances used in the culture solution were: 2 mM L-glutamine, 100 U/mL penicillin, 100 mg/mL streptomycin. Every third-fourth day the cells were fed and reseeded.

Mitoplast Preparation
The experimental methods used for the preparation of fresh mitochondria and subsequent mitoplasts were based on differential centrifugation and hypotonic swelling, respectively, as described in ( [27]). To obtain mitoplasts, the mitochondria were incubated in a hypotonic solution (5 mM HEPES, 100 µM CaCl 2 , pH 7.2) for approximately 1 min, and then a hypertonic solution (750 mM KCl, 30 mM HEPES, and 100 µM CaCl 2 , pH 7.2) was subsequently added to restore the isotonicity of the medium.

Electrophysiological Recordings
The experiments were carried out in mitoplast-attached single-channel inside-out mode for mitoBK channels and also in inside-out mode for the BK channels from the cellular membrane of glioblastoma. The currents were recorded using a pipette of borosilicate glass (Harvard, UK) with a resistance of 7-20 MΩ, which was pulled using the Narishige puller. The glass pipette was filled with an isotonic solution containing 150 mM KCl, 200 µM CaCl 2 , and 10 mM HEPES at pH 7.2 in case of the patch-clamp experiments on the mitoBK channels. In turn, bath and pipette solution contained 130 mM potassium glutamate, 5 mM KCl, 2 mM CaCl 2 , 2.03 mM EGTA, and 8 mM HEPES for the experiments on BK channels from the plasma membrane. In this way, the compositions of the pipette solutions ensured a similar level of Ca 2+ -activation of the investigated BK channels regardless of the differences in sensitivity to calcium ions between both examined channel's isoforms. The current was recorded using a patch-clamp amplifier Axopatch 200B. The currents were low-pass filtered at a corner frequency of 1 kHz and sampled with Clampex software at a frequency of 4.00 kHz (at time intervals of 250.00 µs) in case of mitochondrial patches and 5.00 kHz (at time intervals of 200.00 µs) in case of cell membrane patches. Single-channel currents were recorded at fixed pipette potentials of −60, −40, −20, 20, 40, and 60 mV, at room temperature (22 • C). The measurement error of single-channel currents was ∆I = 1×10 −6 pA (implied by the equipment). Each experimental time series comprised at least N = 2 × 10 5 current values at the applied time resolution of the measurement (the maximal length of recording was N = 6 × 10 6 ). At each value of membrane potential, we recorded time series of single BK channel currents using 3-7 independent patches for each cell membrane type. For each repeating patch-clamp experiment, a fresh patch of mitoplast or cell membrane was used.

Construction of Dwell-Time Series
The construction of dwell-time series corresponding to the subsequent channel states from the original experimental data in a form of time series of single-channel current required choosing appropriate algorithm of event detection i.e., recognition and separation of conducting (functionally open) and non-conducting (functionally closed) state basing on the actual current value. Here we applied the procedure described in [60] to determine the threshold current value (TC) separating the conducting and non-conducting states. The main stages of the TC evaluation are as follows: • plotting the probability density function (PDF) of ionic current approximated by the nonparametric kernel density estimate with Epanechnikov kernel with logarithmic scales on both the horizontal and vertical axes, • finding the intervals, where the power-law scaling is satisfied, • finding the point intersection of the power-law plots (between the unimodal densities representing functionally open and closed states of a channel). The current value corresponding to this point of intersection indicates the TC.

Conditional Mean Dwell-Times of Conducting/Non-Conducting States
We investigated the conditional dwell-time distribution, i.e., the mean opening time duration, provided the previous closing time duration had some particular value (or mean closed time duration provided some duration of the preceding opening event). In this aim, we made use of 2D arrays (similar to histograms), where on the horizontal axis we put the preceding dwell-time duration, and on the vertical axis, we put the average dwell-time duration of a consecutive state. We assumed the time resolution of these arrays of 5 ms.
To understand how to construct such arrays, consider the conditional distribution of opening time duration. Having an opening event, we record the value of the preceding closed time duration. At a resolution of 5 ms, the corresponding array entry for the preceding closed time (the array index calculated as rounding down of t/5 ms) is incremented by the value of the opening duration. In the end, array entries are divided by the number of events accumulated in these entries, to obtain means. The uncertainty is estimated by a traditional uncertainty for the mean of the series aggregated from various realizations of electrophysiological recordings (to account for variability across patches).

Hurst R/S Analysis
The Hurst exponent was calculated on the base of the original Hurst's theory [46]. In brief, the procedure embraces the following stages: first, the original analyzed time series of N elements is divided into n non-overlapped segments of the length s. For each subseries, the mean value E and the standard deviation S are evaluated, and then, the corresponding mean adjusted series and the deviate cumulative series. For each modified subseries (i.e., deviate cumulative subseries) the range R is found as the difference between the maximal and minimal value within a given segment, and it is further rescaled by the standard deviation S. Then, the mean value of the < R/S > is calculated for the current resolution of the division of the original time series (n segments of the length s). Next, the whole procedure is repeated for different lengths s (here, the lengths of subseries are equal to the powers of 2 i.e., s = 2 p , where p = 1, 2, 3, . . . ). Finally, the Hurst exponent is estimated as the slope coefficient of a fitted straight line to a double-logarithmic dependence of form: The values of the Hurst exponent can range from 0 to 1 and can serve as the classification criterion of a given time series in terms of its predictability. If the Hurst exponent takes the value of 0.5, it indicates the purely random behavior of the considered system. In that case, there are no correlations between the subsequent elements of the analyzed series. When 0 < H < 0.5 one can recognize antipersistent characteristics of a given series, which means that the adjacent elements of a given time series tend to switch between the values over and under the mean. If 0.5 < H < 1 that the analyzed time series exhibits trend-reinforcing behavior (its elements are long-term correlated) In this case one can anticipate that a positive increment between the adjacent elements will be followed by a positive one, and a negative increment will be followed by a negative one. The higher values of the Hurst exponent are obtained, the more evident is the trend-reinforcing behavior of the system.

Detrended Fluctuation Analysis (DFA)
The procedure of the DFA analysis [47] (which extends the Hurst approach) determines a level of self-affinity (long-range correlations) of the data even in case of their nonstationarity (statistical characteristics change with time). Here inherent trends are removed at all time scales. First, the signal's profile y i is estimated as the cumulative sum of the series x i with the subtracted average value x Then, the y i is split into n equal non-overlapping segments of size s, which are equal to the powers of two. For each segment, one should find the local trend y p v,i . Here we did it by means of the linear least-squares fit (p = 1). Then, for each v-th segment, the variance is calculated as: and the detrended fluctuation function f (s) for a scale s is given as a square root of the average variance: The detrended fluctuation function is calculated for different lengths of subseries s. Finally, the scaling exponent α is determined from the power-law scaling of the function f (s) with s. It is determined as the slope of the regression line of double-logarithmic dependence log f (s) ∝ α log s.
The α exponent can be interpreted analogously to the Hurst exponent. Namely, if α = 0.5 one can recognize uncorrelated series. If the α lower than 0.5, the analyzed signal is antipersistent. The α values higher than 0.5 characterizes a persistent series. When the α 1 the analyzed signal has 1/f noise characteristics. The values of α above 1 indicate the non-stationarity of the series.
The Hurst R/S and DFA analyses were applied to dwell-time series of subsequent channel states of the length N T = 2 12 (all states) or N Sop/cl = 2 11 (only open or closed states) and to series of single-channel currents of N C = 2 13 samples. The chosen length was dictated by the availability of data (e.g., the length of closed state durations was limited for the series obtained at membrane depolarization). If the original experimental series was sufficiently long, it was divided into non-overlapping subseries of desired length which were analyzed independently.

Multifractal Detrended Fluctuation Analysis (MDFA)
The MFDFA technique can be considered as an extension of the basic DFA technique over the q−th statistical moments of the calculated variance in terms of the scaling function [48] G(q, s) = Consequently, the q-order generalized Hurst exponent H(q) can be calculated from the dependence G(q, s) ∼ s H(q) in a double-logarithmic scale. The mass exponent τ is given by the formula τ(q) = qH(q) − 1. In the following step, the q-order Hölder exponent h(q) is calculated as a derivative of the mass exponent over the argument h(q) = dτ(q) dq . Finally, the q-order multifractal singularity spectrum D(h) (mf-spectrum) is estimated as a Legendre transform of the mass exponent Here, we apply the focus-based formalism of the MFDFA developed by Mukli et al., see Ref. [61] for details. The MFDFA analysis was applied only to the experimental series of single-channel currents of the length N MFDFA = 2 15 . Only in case of this kind of data could we ensure sufficiently long series to detect multiscale multifractal characteristics of the system. The dwell-time series are a composite series and consequently, they have a reduced length which precludes them to be analyzed by means of MFDFA methodology (in the whole range of membrane potentials for a typical duration of patch-clamp recording).

Statistical Analysis
Statistical analysis was performed using Statistica (v13.1PL) and the Pandas [62] statistical package dedicated to Python.
The significance level of α s = 0.05 was selected for the performed statistics. The normality test calculated via Shapiro-Wilk formula does not allow us to confirm the hypothesis about the normal distribution for the majority of analyzed variables. Thus, in further analysis, the non-parametric analog of t-test for the independent samples in the form of Mann-Whitney U statistic was implemented.

Kinetic Features
The samples of experimental patch-clamp traces and the corresponding activation curves (p op (U m )) are depicted in Figure 1.
Our results confirm that the electrophysiological features of different BK channel isoforms stay in a qualitative agreement with each other, both channels are regulated by membrane potential and Ca 2+ concentration. Nevertheless, BK channels from the plasma membrane and their mitochondrial counterparts exhibit different sensitivity to calcium ions quantitatively. Namely, the mitoBK channels are significantly less sensitive to Ca 2+ ions than the cellular ones. The application of the Ca 2+ concentration of 200 µM in case of mitoBK and ca. 18 µM for their cellular analogs imposes a similar level of activation of the calcium sensors for both channel variants and, in consequence, allows to observe a possibly close variability of open state probability at the same range of membrane potential. For both analyzed groups of BK channels, the open state probability increased from below 0.1 to over 0.9 at the applied membrane potentials from −60 mV to 60 mV. As one can see in Figure 1, the applied conditions allowed for a very good concurrence of open state probability at extreme analyzed voltages i.e., −60 mV and +60 mV for both channel variants. Some differences between the voltage-activation curves are visible when comparing the data from mitochondrial and cellular BK channels in glioblastoma at medium values of voltage. Namely, the BK channels from mitochondrial membranes are more sensitive to membrane depolarization than their analogs from the plasma membrane. We think that these discrepancies can stem from the different interactions between the voltage and/or calcium sensors and the channel's pore-gate domain in case of the cell membrane and mitochondrial BK channel variants, or differences in lipid and protein composition of membranes which can further alter the channel gating. In fact, these differences in BK and mitoBK channel functioning were further investigated by means of several techniques of correlation analysis as a merit of this work. What is here also important, the applied measures of short-and long-term correlations do not directly depend on the average conductance of the channel. They refer to the structure of the signal and its complexity.

Conditional Mean Dwell-Times of Conducting/Non-Conducting States
The results of the conditional dwell-time duration distributions differ for the mitochondrial and cellular channels (Figures 2 and 3).

As one can see in Figures 2 and 3 at negative voltages long dwell-times of open states tend to be next to relatively short closed dwell-times, and vice versa. This kind of behavior is a manifestation of short-range correlations between the subsequent durations of consecutive channel states.
To understand (without going into too many details) the meaning of the above analysis, let us consider some simple gating models. What is here well worth mentioning, the simplest possible description of the channel kinetics by a two-state Markovian model (see Figure 4a) is not able to describe any correlations between the dwell-times of channel states. The correlations can arise when channel kinetics is modeled by at least two substrates within the manifolds of open and closed channel states, and the open and closed states of the different average lifetime are connected to each other [63,64]. This can be manifested in a form of loops within the kinetic scheme which is shown in Figure 4b presenting a possibly simple representation of BK channel gating by 3 open and 5 closed states at fixed experimental conditions (that kind of modeling can be further extended to multiple state Markovian models to describe broad ranges of voltage and Ca 2+ concentrations that affect gating). In the general case, 3-4 open states and 5-6 closed states can model BK channel gating [65]. There are also several common multistate extensions of the Markov model [66,67], that allow for description of the gating kinetics in a broad range of voltages and concentration of calcium ions. In this paper, we do not evaluate the number of states or its connectivity based on the experimental results, but we would like to make the following remarks. Firstly, the kinetics of the channels from both mitochondrial and cellular membranes should be described by a relatively complex model with a large net of the open and closed substrates. Secondly, the recognized differences in conditional dwell-times of open and closed states are a direct consequence of the changes in the Markov diagram describing the considered channel isoforms. One can expect different sizes of the open and closed state manifolds or different rate constants for state transitions. In turn, the differences in Markov state diagrams are consequences of structural differences of the channel variants expressed in cell membranes and mitochondria that may affect interactions between pore-gate and the voltage/calcium sensors and gating dynamics (also different additional factors influencing conformational diffusion space of the BK channel protein may occur in mitochondrial and cellular membranes e.g., another lipid and protein partners which can interact with the BK channel).

Hurst R/S and DFA Analyses
Due to the fact that the DFA method can be considered as a generalization of the R/S methodology, we decide to provide the results obtained by both methods in one section. Admittedly, the DFA reveals more evidently the long-range memory effect than the R/S analysis. Still, the general inference, based on the Hurst and scaling exponents describing the single-channel activities from different patches is that both single-channel currents, as well as dwell-time series, unravel long-range memory effect for the BK channels' activity regardless of the membrane type (mitochondrial or cellular), applied voltage, or whether we analyze the total signal or its components corresponding to functionally open or closed states. This makes that the existence of long-range correlations an inherent feature of gating dynamics.
It is difficult to indicate unequivocally the origins of the long-range memory in the channel's activity. The discrete several-states Markovian models (3-, 4-and 11-state models) are not able to reproduce long-term persistence measured by the Hurst exponent [49,50]. Looking for an explanation of the causative factors inclining the Hurst memory effect, we support the previously formulated hypothesis [54], where the possible source of long-term correlation is related to the presence of a component side-process occurring at a larger time scale in relation to the conformational diffusion of a channel protein but still deeply influencing the diffusive space of the channel gate (like membrane fluctuations).
The recognized changes of the standard and generalized Hurst exponents with the increasing U m are relatively weak and non-monotonic. This may be surprising considering that at deep membrane hyperpolarization short openings tend to follow long-lasting closed states. At medium values of membrane potential both dwell-time durations equilibrate. Finally, as the voltage rises up to a high depolarization the dwell-times of channel states reach an extremum, where short closings tend to follow long-lasting open states. However, as stated above, it's not the Markovian part of the channel dynamics that drives the long-term correlations.
From the reason that there is no clear tendency between the Hurst exponents and U m we grouped the sets of H (R/S) and α (DFA) exponents obtained for a given membrane-type at all the applied voltages and then compared medians of H and α representing mitochondrial and cellular patches in three variants: total signal, and the series of open and closed states separately.

Analysis of Dwell-Time Series of Channel States
The analysis of the dwell-time series of subsequent channel states by the Hurst R/S and DFA analyses allows for clear discrimination between the mitochondrial and cellular channels ( Figure 5).
The memory effect is evidently lower for the channels from mitochondrial patches than for the cellular ones ( Figure 5 Except for one case of the mean Hurst exponent values calculated by the R/S method for the closed state (p = 0.645 for Mann-Whitney U test, see Figure 5f) all the differences were statistically significant (p < 0.05) when one compares the dynamics of channels from mitochondrial and plasma membranes.
Thus, one can observe that application of the R/S and DFA methods gave almost consistent results, but the recognized relation between the H mito vs. H cell is slightly better pronounced in the case of the DFA methodology (which possibly fits better to the structure of analyzed signals where the underlying statistics and dynamics can be non-stationary). The larger long-term memory effect between the dwell-times of subsequent channel states in the case of channels from cellular membrane indicate that in this case, the mechanism of switching between the states of different lifespan favors more evidently the pattern that the long-lasting sojourns in subsequent states should be followed by a long-lasting ones, and the short-lasting channel states should be followed by states of short durations. What is interesting this feature of the system is visible as a general rule occurring over a range of time scales (as dictated by the R/S or DFA methodology) except the short time scales, like in the case of the calculated conditional open and closed dwell-times (Figures 2 and 3), where it disappears for both the mitoBK channels and their analogs form plasma membrane (this feature seems to be lost in the "Markovian" noise of the channel). The relation between the Hurst exponents, i.e., H mito < H cell may be interpreted according to the aforementioned framework of membrane thickness fluctuations [54].
Namely, a part of the fluctuations' range may be insufficient to alter the diffusive space of the channel gate. In the case of the inner mitochondrial membrane, this fraction could be greater than in the case of the plasma membrane. In consequence, the fluctuations of the mitochondrial membrane seem to affect the mitoBK gating machinery to a smaller extent in comparison with the effects of plasma membrane fluctuations on the cellular variants of the BK channels. What is also interesting, the deviations of the Hurst exponents describing BK channels from the cell membrane are significantly higher than those characterizing mitoBK channels ( Figure 5). Thus, the memory effect has a significantly higher versatility for the plasma membrane channels-they are able both to present a strongly long-range correlated pattern of gating as well as antipersistent one. The presence of these outstanding subgroups of channel behaviors implicates that the population of cellular BK channel is not as uniform in terms of the characteristics of ionic signaling as the mitochondrial one.

Analysis of Single Channel Currents
The persistent behavior of the system is even more evident (higher H and α exponents) in case of the analysis of single-channel current traces than the dwell-time series describing lifespan of subsequent channel states (analogously, like in the [55]), as shown in Figures 5 and 6. Considering the results, the statistically significant differences of α exponents were obtained for all analyzed currents (p < 0.05) but only in the case of the DFA method (Figure 6a-c), in contrast to the R/S where there is a lack of sufficient grounds for rejecting the null hypothesis of no statistical differences between the compared variables: p = 0.440 in the case of the total signal, p = 0.239 for the currents recognized as channel openings and p = 0.123 for the fluctuations of recorded current during the closed state of a channel ( Figure 6 d-f). These results can confirm the hypothesis that in relation to the R/S technique, DFA is a method that allows one to determine more subtle differences between the series characterized by a high degree of complexity. For the currents corresponding to the closed channel's states we obtained H mito = 0.72, α mito = 0.94 and H cell = 0.73, α cell = 0.81, and for the ones corresponding to the conducting ones: H mito = 0.68, α mito = 0.86 and H cell = 0.67, α cell = 0.73, where the exponents are given as medians from the data obtained at all membrane potentials.
As one can see the switching between the channel states of different pore conductance is realized in a more regular way (persistent) in the case of mitochondrial channels than for the ones located in the cell membrane. This effect can have its representation in the conditional dwell-time curves (Figures 2 and 3) where in most cases the dwell-times of channel sojourns in the subsequent states reach higher values, so the mitochondrial channel tends to occupy conformations of similar conductance classified as one macrostate (conducting i.e., open or non-conducting i.e., closed) for a long time.
Similarly to the previous case when the dwell-times were analyzed, investigations of the long-range correlations within the single-channel currents also indicate that the total number of Hurst exponent outliers and the corresponding deviations are higher in the case of the cellular BK channels than the mitochondrial ones ( Figure 6), which confirms the higher variability of gating behavior among cellular channels (a possibly more complex mixture of channel isoforms present in cell membrane than in the mitochondrial membrane; however the channels of outstanding characteristics are relatively rare-they do not affect the median considerably).

MFDFA Analysis of Experimental Single-Channel Recordings
The final correlation analysis performed in this work was the MFDFA analysis of single channel currents. This methodology allowed us to indicate that the series of channel currents have multifractal characteristics and are nonrandom but caused by the orderly process exhibiting long-range correlation features regardless of the applied voltage or membrane type (cellular or mitochondrial). As one can see in Figure 7, the obtained multifractal spectra corresponding to BK channels from the inner mitochondrial membrane are well separated from the ones describing the activity of channels from plasma membranes. Only at strong membrane potential (U m = ±60 mV) the gating dynamics of mitochondrial and cell-membrane BK channels seem to be less diversified. Figure 7. Average multifractal spectra obtained for BK channel from mitochondrial and cellular patches at different membrane potentials. The blue spectra characterize the mitoBK channel's activity; the orange spectra are assigned to the BK channels from the cell membrane. The upper panels (a-c) present the results obtained at hyperpolarization; the lower graphs (d-f) correspond to the channels' characteristics at membrane depolarization.
Due to the fact that the spectral width (∆) has the most influential interpretative potential than the other parameters of the spectra, we restricted the presentation of the obtained results to the ∆, as shown in Figure 8. The obtained results show that the spectral width is evidently higher in the case of the BK channels from mitochondrial patches than the ones from the cellular membrane, which suggests a higher complexity and entropy of the activity of the mitoBK channel [68,69]. All the presented differences between the mf-spectrum width of cellular and mitochondrial series reached the statistical significance for the Mann-Whitney U test (p < 0.05).
According to Boltzmann's definition of entropy (where the entropy is proportional to the natural logarithm from the number of real microstates of a considered system), these results ought to indicate a larger number of the attainable substates (stable conformations) in the case of the mitochondrial channel than its cellular analog. A hypothetical biological interpretation of this phenomenon is that the mitochondrial-membrane channels exhibit higher affordability to switch between the complex manifold of channel's conformations, which can potentially improve their adjustment to an ever-changing environment.
The larger width of a spectrum at both strong membrane depolarization and hyperpolarization comparing with the ones obtained at moderate membrane potentials in both cases which can be a general feature of the BK channel dynamics regardless of its splice variant (cellular or mitochondrial), see Figure 8. From the definition of Boltzmann's entropy, it seems that that recognized greater complexity of the signal at highly positive and negative potentials in comparison with the data obtained at membrane potentials close to zero should stem from a possible increase in the number of attainable channel substates with the absolute value of voltage. The qualitative agreement between characteristics of mitochondrial channels and the ones from the cell membrane suggests that this feature should stem from the interactions between the pore-gate domain and the voltage sensor (the structure of the membrane-spanning domains should be identical in both BK channel isoforms). In turn, the quantitative differences can originate mainly from different lipid surroundings within the plasma and mitochondrial membrane which can affect the BK channel conformational space.

Analysis of Randomized Series
To confirm whether the long-range memory effect recognized for the experimental data is its intrinsic feature, the investigated series were shuffled, and then the Hurst and scaling exponents as well as the MFDFA spectra were calculated for the randomized data sets. The randomization of analyzed series describing single-channel activity-both currents and dwell-times of channel states, caused the vanishing of the long-range memory effect measured by the Hurst exponent (H ≈ 0.50 ± 0.01 for all investigated cases). In turn, during the DFA analysis of shuffled series we encountered quite an interesting feature of the data. The randomization procedure applied to the analyzed series of single channel currents caused vanishing of the long-range correlations (α ≈ 0.53 ± 0.02 for all investigated patches). Nevertheless, a quite surprising result was obtained by shuffling the series of dwell-times of subsequent channel states. Namely, it does not fully erase the long-range correlations, but it allows us to significantly reduce them (all α values decreased by about 0.10). The possible explanation of this effect is that the distribution of dwell-times are not normal, and the dwell-times itself are conjunctive characteristics of the experimental signal. It turns out that the current work is not a sole report in which the long-range memory effect does not fully vanish after the shuffling procedure. In [49,50], as here, the Hurst exponent for the shuffled series was not exactly 0.50 (but about 0.57 ± 0.02). Still, most importantly, the application of statistical tests indicated that the Hurst exponent for shuffled data was significantly different from the mean value of the original experimental Hurst coefficient.
Shuffling operation of the series analyzed by the MFDFA method resulted in considerable changes of the spectral distributions and consequently affected significantly the spectral parameters (the spectra of shuffled data are about twice narrower). Reduced multifractality after mixing of the data gives clear evidence that the recognized fractal nature of the examined signals is its inherent feature and the system dynamics can be considered as an orderly process exhibiting long-range correlations.

Conclusions
In this work, the existence of dynamic diversity within the BK channels from the plasma membrane and the inner mitochondrial membrane in human glioblastoma cells is examined.
Our results indicate that the mitoBK channels and their cell membrane analogs exhibit different sensitivity to calcium ions, which was observed or postulated by the other authors [18,29,[43][44][45].
Probably the weakened sensitivity of mitochondrial BK channels stems from the structure-function relationship of the particular splice variants present in mitochondria or specific protein-protein and protein-lipid interactions within the inner mitochondrial membrane, which can deeply affect channel functioning and allows for the adaptation of channel transport efficacy to the local conditions; in particular relatively high Ca 2+ concentration in mitochondria (which act as calcium reservoirs). It can be anticipated that the ligand sensors in both cases (mitochondrial and cellular BK channels) can exhibit other dynamics of calcium-binding. Probably, the Ca 2+ -sensors from the BK channels from plasma membrane can bind the calcium ions for a longer time than their mitochondrial analogs. This allows the BK channels from the plasma membrane to reach a low Ca 2+ concentration (where the waiting time for the Ca 2+ ions to approach the Ca 2+ -free sensors is relatively large) a similar level of Ca 2+ channel activation like the mitoBK channels at large Ca 2+ concentration.
In the experimental part of our research, we attempted to adjust the external conditions to impose possibly close levels of voltage and calcium activation of the analyzed mitochondrial and cell-membrane channels to enable detection of possible differences in channel gating between these two isomeric BK channel variants.
We laid a special emphasis on the description of short-and long-range correlations within the experimental data describing single-channel activity. The obtained results allowed us to indicate both generic features of the gBK channels which are common for their mitochondrial and cellular variants (they are based on the qualitative agreement between the evaluated characteristics) and also point out the differences that allow us to discern both groups of channels (quantitative differences of evaluated measures between mito-and cellular channels). Among the inherent features of BK channel activity referring to the correlations in system dynamics are: • short-range correlations between the lifespan of subsequent channel states confirmed by the analysis of the conditional dwell-time duration as shown in Figures 2 and 3, which suggest the connectivity between the channel's substates of different longevity (i.e., the number of routes between the substates within the manifold of functionally open or closed states) [63,64], • long-range persistence indicated by Hurst exponents calculated by the standard R/S and DFA algorithms for single-channel currents ( Figure 6) and the series of dwell-times of subsequent channel states ( Figure 5), • multifractal complex dynamics indicated by the MFDFA analysis (Figure 6), where the spectral width exhibits an analogous dependence on the membrane potential for both types of BK channels variants (broader spectrum at highly hyper-or depolarized membranes than at the U m close to zero) (Figure 8), suggesting similar changes of the conformational diffusive space with U m .
These common characteristics can originate from the interactions and functioning of transmembrane segments of BK channels which should be almost identical for mitochondrial and cellular splice variants of BK channels in glioblastoma.
The performed analyses allowed also to indicate some of the details that distinguish the dynamics of mitochondrial and cell-membrane BK channels variants. The anticipated presence of some structural differences between the analyzed groups of BK channels resulted in e.g., the visible different Ca 2+ sensitivity. For this reason, the experimental data used in this work were obtained at such conditions which ought to effectively attenuate these differences and impose a similar level of calcium activation of the channel. In consequence, we could observe the effects of interactions between the channel's sensors within the channel isoforms from mitochondrial and cell membranes, and detect the possible subtle differences in the system's dynamics.
The first evidence of such dynamical diversity was observed during the investigations of the conditional dwell-times of open and closed states of the channels. As mentioned in Results and Discussion Section this is a method for probing the Markov state diagram of the channel. Different states display different rate constants, so after estimating the mean time duration following a particular dwell-time, one may anticipate the connectivity between two states of the diagram. For example if one considers an open state duration that follows a closed state that lasts for 10 ms, one should investigate only the opening state, that is accessible from a closed state of liveliness 10 ms, but not from a state of liveliness of e.g., 30 ms. Obviously, if there are structural or functional differences in the channel, we expect the differences in the Markov state model that would describe such a channel. We have found differences in the conditional dwell-time distributions that indicate changes in the underlying Markov model, so the existence of the structural or functional differences between mitochondrial and cellular channels appears to be indirectly confirmed in that way. The questions that remain are: what is the origin of these differences? Is it only genetic splicing or folding conditions within the same genetic sequence? Possibly, the hypothetical slight changes of the dynamics of the Ca 2+ binding, and further allosteric interactions between the sensors and pore-gate domains can also modulate gating dynamics (and the examined correlations). Another process that can influence gating kinetics can be connected with different compositions of cellular and mitochondrial membrane, and consequently other lipid-protein and protein-protein interactions [21][22][23][24][25]. An important side-effect of the conditional dwell-time analysis is that it should give good constraints to estimate the rate constants of the many-state BK channel models. This is, however, one of the problems that we would like to address in further research.
The analysis of long-range correlations indicates a larger memory effect within the series of dwell-times of subsequent channel states in the case of the BK channels from the cell membrane than the mitochondrial ones. This relation is observed irrespective of whether all channel states are considered or the openings and closings separately. Opposite relations are obtained when the series of channel currents are analyzed. The single-channel currents and the dwell-time series are however separate characteristics of the system, so such inversion of the relations is acceptable. What is worth mentioning it is that for all analyzed isoforms of the channel the long-range correlated characteristics of the non-conducting states of a channel determined the properties of the series describing the total trace of channel's activity (higher H and α values at channel's closings than openings in most cases), which is in good agreement with the former reports [55,56,59].
Summing up, the obtained results indicate that the splicing sequence or interactions with other membrane components corresponding to a particular location can define the conformational machinery of channel gating. This in turn allows us to adjust the channel functioning to a local environment and enable effective fulfilling of the physiological tasks of the channel in a given location. The recognized differences in gating dynamics between the channel variants can be represented as different connectivity, rate constants of transitions between the channels' substates and/or even the number of attainable substates in appropriate kinetic schemes describing channel-gating (as in a popular Markov-type modeling), as postulated in this work. Considering the significance and potential utility of the obtained results this work shows that the methods of correlation analysis can serve as a useful tool to discern sets of experimental data describing the activity of channels that are structurally and functionally similar, as different groups of isoforms of the same channel type. Thus, the measures of short-and long-range correlations can be a basis for some classification algorithms in artificial intelligence approaches aimed to determine patterns within the experimental data describing ion-channel activity with or without the use of specific channel modulators. Their implementation in AI techniques (along with the measures of gating kinetics) can simplify the choice of the specific activating/inhibiting substances which modulate the channel functioning in the desired manner. The results obtained here can provide a theoretical background for further investigations on the potential modulators of the mitochondrial BK channels as well as their cellular analogs. As our outcomes suggest differences in gating machinery between the mitoBK and the BK channels from the plasma membrane in glioblastoma cells they allow us to infer that there should be a possibility of a highly specific modulation of only one group of them. Investigations of such kind of modulators could be one component in the development of new treatment strategies against gliomas which constitutes an important venture in future research.
Thinking beyond the glioblastoma treatment, one could pose also a general question about the possible utility of active substances which could affect the memory effects within the channel gating. Admittedly, the biological meaning of long-range correlations is still not fully understood. However, according to the definitions of the measures calculated in this work, this kind of therapeutics could play an important role in the regulation of biological processes whose characteristics reach a long-range in time. A reasonable candidate of diseases meeting this criterion are arrhythmia or epilepsy, for example. The potential meaning of BK channels in these diseases was summarized in the works [70,71] however the role of channel gating memory in their pathogenesis has never been discussed yet. One can anticipate that these diseases can be correlated with the abnormalities in long-range correlations in channel gating, which has to be confirmed in further studies. Funding: Agata Wawrzkiewicz-Jałowiecka was supported by the Silesian University of Technology Grant for young researchers BKM-528/RCh4/2020 (statute project).