Computational prediction of drug response in short QT syndrome type 1 based on measurements of compound effect in stem cell-derived cardiomyocytes

Short QT (SQT) syndrome is a genetic cardiac disorder characterized by an abbreviated QT interval of the patient’s electrocardiogram. The syndrome is associated with increased risk of arrhythmia and sudden cardiac death and can arise from a number of ion channel mutations. Cardiomyocytes derived from induced pluripotent stem cells generated from SQT patients (SQT hiPSC-CMs) provide promising platforms for testing pharmacological treatments directly in human cardiac cells exhibiting mutations specific for the syndrome. However, a difficulty is posed by the relative immaturity of hiPSC-CMs, with the possibility that drug effects observed in SQT hiPSC-CMs could be very different from the corresponding drug effect in vivo. In this paper, we apply a multistep computational procedure for translating measured drug effects from these cells to human QT response. This process first detects drug effects on individual ion channels based on measurements of SQT hiPSC-CMs and then uses these results to estimate the drug effects on ventricular action potentials and QT intervals of adult SQT patients. We find that the procedure is able to identify IC50 values in line with measured values for the four drugs quinidine, ivabradine, ajmaline and mexiletine. In addition, the predicted effect of quinidine on the adult QT interval is in good agreement with measured effects of quinidine for adult patients. Consequently, the computational procedure appears to be a useful tool for helping predicting adult drug responses from pure in vitro measurements of patient derived cell lines. Author summary A number of cardiac disorders originate from genetic mutations affecting the function of ion channels populating the membrane of cardiomyocytes. One example is short QT syndrome, associated with increased risk of arrhythmias and sudden cardiac death. Cardiomyocytes derived from human induced pluripotent stem cells (hiPSC-CMs) provide a promising platform for testing potential pharmacological treatments for such disorders, as human cardiomyocytes exhibiting specific mutations can be generated and exposed to drugs in vitro. However, the electrophysiological properties of hiPSC-CMs differ significantly from those of adult native cardiomyocytes. Therefore, drug effects observed for hiPSC-CMs could possibly be very different from corresponding drug effects for adult cells in vivo. In this study, we apply a computational framework for translating drug effects observed for hiPSC-CMs derived from a short QT patient to drug effects for adult short QT cardiomyocytes. For one of the considered drugs, the effect on adult QT intervals has been measured and these measurements turn out to be in good agreement with the response estimated by the computational procedure. Thus, the computational framework shows promise for being a useful tool for predicting adult drug responses from measurements of hiPSC-CMs, allowing earlier identification of compounds to accurately treat cardiac diseases.

the computational procedure appears to be a useful tool for helping predicting adult drug responses from pure in vitro measurements of patient derived cell lines.

Author summary
A number of cardiac disorders originate from genetic mutations affecting the function of ion channels populating the membrane of cardiomyocytes. One example is short QT syndrome, associated with increased risk of arrhythmias and sudden cardiac death. Cardiomyocytes derived from human induced pluripotent stem cells (hiPSC-CMs) provide a promising platform for testing potential pharmacological treatments for such disorders, as human cardiomyocytes exhibiting specific mutations can be generated and exposed to drugs in vitro. However, the electrophysiological properties of hiPSC-CMs differ significantly from those of adult native cardiomyocytes. Therefore, drug effects observed for hiPSC-CMs could possibly be very different from corresponding drug effects for adult cells in vivo. In this study, we apply a computational framework for translating drug effects observed for hiPSC-CMs derived from a short QT patient to drug effects for adult short QT cardiomyocytes. For one of the considered drugs, the effect on adult QT intervals has been measured and these measurements turn out to be in good agreement with the response estimated by the computational procedure. Thus, the computational framework shows promise for being a useful tool for predicting adult drug responses from measurements of hiPSC-CMs, allowing earlier identification of compounds to accurately treat cardiac diseases.
1 Introduction 1 Short QT (SQT) syndrome is a cardiac channelopathy characterized by an 2 abnormally short duration of the QT interval of the patient's electrocardio-3 gram (ECG) [1,2,3]. The syndrome is associated with increased risk of 4 atrial fibrillation, ventricular arrhythmias and sudden cardiac death [4,5] 5 and was first described by Gussak et al. in 2000 [1]. The first identified SQT 6 subtype, termed SQT1, results from an increase in the transmembrane potas-7 sium current, I Kr , caused by a mutation in the I Kr encoding gene KCNH2 8 [6]. Later, several additional subtypes of SQT syndrome have been identi-9 fied, originating from mutations in other genes, including gain of function 10 alterations in genes encoding the potassium channels responsible for the I Ks 11 and I K1 currents and loss of function alterations in genes encoding calcium 12 channels [2,3]. 13 Because of the lethal risks associated with the syndrome, there is an ur-14 gent need for effective pharmacological therapies. One modern approach is 15 to identify compounds that can make the electrophysiological properties of 16 cardiomyocytes (CMs) populated with SQT mutated channels more simi-17 lar to the electrophysiological properties of CMs populated with wild type 18 (WT), unmutated, channels (see e.g., [7,8,9,10,11]). For SQT1, which is 19 characterized by an increased I Kr current, drugs inhibiting the I Kr current 20 have been proposed as possible candidates (see e.g., [11,12]). However, the 21 effect of a drug on WT I Kr channels can be very different from the effect 22 on mutated channels. For example, the I Kr blockers sotalol and ibutilide 23 have been shown to be ineffective for SQT1 patients [12]. This may be ex- 24 plained by the fact that these drugs primarily affect the inactivated state 25 of the I Kr channels, and the SQT1 mutation impairs the inactivation of the 26 channels [6, 13,14]. On the other hand, quinidine, which affects both the 27 open and inactivated states of I Kr channels has proven to be more effective 28 for SQT1 patients [15]. These examples demonstrate that investigations into 29 drug effects on specific SQT mutations are needed in order to find appropriate 30 pharmacological treatments. 31 A promising prospect in the quest for suitable drugs for SQT syndrome 32 is the development of human induced pluripotent stem cells (hiPSCs) (see 33 e.g., [16,17,18,19]). These cells can be generated from individual patients 34 and differentiated into a large number of different cell types, including CMs. 35 Thus, cells can be generated from SQT patients, allowing for investigations 36 of drug effects for CMs populated by ion channels affected by the patient 37 specific mutations. Indeed, in [9,11], several drugs have been tested on 38 hiPSC-derived CMs (hiPSC-CMs) from an SQT1 patient, revealing possible 39 promising drugs. 40 However, a difficulty with using hiPSC-CMs in drug testing applications 41 is that the electrophysiological properties of hiPSC-CMs differ significantly 42 from the electrophysiological properties of adult CMs (see e.g., [20,21]). In 43 general, hiPSC-CMs are recognized as electrophysiologically immature, with 44 properties more similar to those of fetal CMs. These differences imply that 45 the drug response observed for hiPSC-CMs might not be the same as the 46 corresponding drug response for adult native CMs. 47 Mathematical modeling has been proposed as a possible tool to help trans- 48 late the drug response of hiPSC-CMs to the drug response of adult CMs (see 49 e.g., [22,23,24,25]). The development of mathematical models of the dynam-50 ics underlying human cardiac action potentials is an active field of research, 51 and a large number of models have been developed, including models for 52 adult undiseased ventricular CMs [26,27], atrial CMs [28,29], hiPSC-CMs 53 [30,31] and CMs affected by mutations [32,33]. In these models, changes in 54 the membrane potential are represented by individual transmembrane ionic 55 currents (see e.g., (1) below), and the models can therefore be useful for 56 investigating how drug effects on individual ion channels affect the compos-57 ite dynamics of the full action potential. Furthermore, the action potential 58 models can be combined with spatial models of electric conduction (see e.g., 59 [34,35,36,37,38]) to provide insight into drug effects on cardiac conduction 60 properties and mechanisms of arrhythmia. 61 In previous studies, we have applied a procedure based on mathematical 62 action potential modeling to estimate drug effects on individual ion channels 63 from measurements of hiPSC-CMs and predict adult CM drug responses 64 [22,23]. This procedure is based on the assumption that the function of 65 individual proteins, e.g. ion channels, is the same for hiPSC-CMs and adult 66 CMs, and that only the density of the proteins differs between hiPSC-CMs 67 and adult CMs (see Figure 2). From this assumption, it follows that the 68 effect of a drug on an individual ion channel is the same for hiPSC-CMs and 69 adult CMs. Assuming that we have correctly identified the drug effect on 70 individual channels in the hiPSC-CM case, we can therefore directly translate 71 this effect to the adult case by inserting the inferred mechanisms into a model 72 for adult cells. 73 The aim of the present study is to use this computational procedure to 74 predict drug effects for adult SQT1 CMs based on measurements of drug ef-75 fects on the action potential of SQT1 hiPSC-CMs from [9,11] and to extend 76 these results into prediction of patient QT changes. Our overall computa-77 tional pipeline is depicted in Figure 1. We consider the four drugs quinidine, 78 ajmaline, mexilietine and ivabradine, shown to potentially be useful for SQT1 79 patients in [9,11] and show predicted drug responses for the drugs on the 80 ventricular action potential and QT interval of adult patients. We validate 81 our pipeline using data for quinidine, where measurements of the drug effect 82 on the QT interval have been conducted for adult patients [15]. The pre-83 dicted drug response turns out to be in good agreement with the measured 84 drug effect, indicating that the computational procedure could be useful for 85 predicting adult drug responses from measurements of hiPSC-CMs. 86

87
In this section, we describe the methods applied in this study. We start by 88 describing the basic modeling assumptions underlying the approaches used 89 for computational identification of drug effects and mapping of drug effects 90 from hiPSC-CMs to adult CMs. Then, we describe details of the applied 91 action potential model and inversion procedure. Finally, we describe the 92 Figure 1: Illustration of the computational pipeline. 1) Biomarkers from the cardiac AP are taken from hiPSC-CMs under drug testing. 2) These biomarkers from dose escalation studies are inverted into an SQT1 model of the AP of hiPSC-CMs. Inversion into a matched model provides determination of drug effects on specific channels [23]. 3) The drug effects determined in 2 are inserted into a model of adult CMs with the same SQT1 mutation to give a prediction of drug effect in mature CMs [22]. 4) The adult CM model is converted into pseudo QT waveforms for prediction of QT segment changes in SQT1 patients under the estimated effect of the drug. approach used to estimate the QT interval in the adult case. Note that 93 the majority of these methods are to a large extent based on the methods 94 described in [23]. 95 2.1 Action potential models 96 In this section, we describe the framework for modeling the action potentials 97 of hiPSC-CMs and adult CMs, with and without the SQT1 mutation, and the 98 relationship between these models. Three of the main modeling assumptions 99 are illustrated in Figure 2, building on the approaches of [22,23]. In the 100 next subsections we will explain these assumptions and demonstrate how 101 they affect the action potential modeling. In the action potential model, the membrane potential, v (in mV), is governed 104 by an equation of the form where I j are different membrane current densities. These current densities 106 can be expressed on the form where N j is the number of proteins of type j on the membrane, A is the 108 area of the membrane, C m is the specific membrane capacitance, and i j is 109 the average current through a single protein of type j. For currents through 110 voltage-gated ion channels, this i j is typically given on the form where g 0,j is the conductance through a single open channel, o j is the open 112 probability of the channel, and E j is the equilibrium potential of the channel. 113 In this case, it is common to introduce combined parameters of the form and write (2) on the form  Figure 2: Illustration of the assumptions underlying the computational maturation approach. 1) The density of different types of ion channels (and other membrane or intracellular proteins) may differ between hiPSC-CMs and adult CMs, but the function of the individual channels is the same. In the model, the density difference is represented by the parameter λ.
2) The SQT1 mutation affects the individual I Kr channels in exactly the same manner for hiPSC-CMs and adult CMs. In the model, the mutation is represented by an adjusted model for the open probability, o Kr .
3) The effect of a drug on a single protein is the same for hiPSC-CMs and adult CMs. In the model, the drug effect is represented by the parameter ε. A number of electrophysiological properties have been shown to differ be-118 tween hiPSC-CMs and adult CMs (see e.g., [20,21]). In addition, the elec-119 trophysiological properties of different samples of WT hiPSC-CMs often vary 120 significantly (see e.g., [39,40]). We assume that these differences in electro-121 physiological properties (both between samples of hiPSC-CMs and between 122 adult CMs and hiPSC-CMs) are due to:  However, we assume that the function of the individual membrane and in-129 tracellular proteins are the same for the different WT cases.

130
Considering the model (1)-(2), these assumptions imply that the density 131 ρ j = N j A may differ between cells, but that i j remains the same. Parameter-132 izing a model to specific set of measurements can then be accomplished by 133 adjustments of the densities represented by adjustment factors λ j such that 134 where ρ * j is the default density of proteins of type j. Incorporating these adjustment factors into the model (1), specific adult and hiPSC-CM versions of the model are given by respectively, where λ A j and λ hiPSC j specify the protein densities in the two 135 cases. For currents through ion channels, these adjustment factors can be 136 incorporated by adjusting the conductances, g j (see (4)-(5)). Similar adjust-137 ment factors can be set up for the density of intracellular proteins and the 138 cell geometry, see [23]. A list of the maturation-dependent parameter values 139 of the base model are given in Table S4  We assume that a mutation affects an individual ion channel in exactly the 143 same manner for hiPSC-CMs and adult CMs. SQT1 is known to affect the 144 I Kr current [41,3], and we assume that the effect on the current can be rep-145 resented by adjusting the model for the open probability, o Kr (see (3)). More 146 specifically, we adjust the voltage dependence of the steady state inactivation 147 gate, x Kr2 by shifting the inactivation towards more positive potentials, as 148 described for the SQT1 mutation N588K in [3]. Because we assume that the 149 mutation affects an individual I Kr channel in the same manner for hiPSC-150 CMs and adult CMs, we apply the same adjustment of the model for o Kr in 151 the hiPSC and adult SQT1 cases. The third assumption illustrated in Figure 2 is that since the function of a 154 single ion channel is the same for hiPSC-CMs and adult CMs, the effect of a 155 drug on an individual ion channel will also be identical for hiPSC-CMs and 156 adult CMs. This means that if we are able to determine the effect of a drug 157 on an individual channel in the hiPSC-CM case, we can find the effect of the 158 same drug in the adult case by incorporating the same single channel drug 159 effect in the adult version of the model.

160
Following [23], we use a simple IC 50 -based modeling approach to represent 161 the effect of a drug. That is, we assume that the average single channel 162 current through a channel j in the presence of the drug dose D is given by 163 Here, i j (0) is the average single channel current when no drug is present and 164 ε j (in µM −1 ) represents the effect of the drug, defined as where IC j 50 (in µM) is the drug concentration that blocks the current j by 50%. Incorporating drug effects into the hiPSC-CM and adult CM models where we note that ε j is the same in the two versions of the model.

166
In this study, we will use this modeling framework to estimate the drug 167 effect (in the form of ε values) of drugs based on action potential measure-168 ments of hiPSC-CMs with the SQT1 mutation N588K. Next, we will insert 169 the estimated ε values into a model for adult ventricular cells with the same 170 mutation to estimate the effect of the drug for an adult patient. In order to represent the action potentials of adult ventricular CMs and 173 hiPSC-CMs (both with and without the SQT1 mutation N588K), we use 174 a slightly modified version of the base model formulation from [23], follow-175 ing the form (1)- (2). Three main adjustments have been incorporated into 176 the current version of the model. First, we have modified the steady-state 177 values of the gating variables of the I Kr current to better fit measured I Kr 178 currents in the WT and SQT1 cases from [14] (see Figure 4). Second, we 179 have extended the temperature-dependence of the model by including Q 10 180 values for a number of the currents, following [27]. This was done because 181 we consider two different temperatures in our computations. In the adult 182 case, we assume body temperature (310 K), while for the hiPSC-CM case, 183 we assume room temperature (296 K) since the considered hiPSC-CM data 184 (from [9,11]) are obtained at room temperature. Third, we have incorpo-185 rated a dynamic model for the intracellular Na + concentration. This was 186 done in order to make the frequency dependent QT interval changes more 187 physiologically realistic (see Figure 10). The full base model formulation is 188 found in the Supplementary Information. The system of ordinary differential 189 equations (ODEs) defined by the base model is solved using the ode15s solver 190 in Matlab. In the adult case, we use 1 Hz pacing unless otherwise specified. Furthermore, 193 we run each simulation for at least 500 pacing cycles after each parameter 194 change in order to obtain new stable solutions before recording the action 195 potentials (and pseudo-ECGs, see Section 2.4 below).

196
In the hiPSC-CMs case, we use 0.2 Hz pacing, as specified in [11]. As a 197 compromise between the need to obtain new stable solutions after each pa-198 rameter change and the need for reducing the computing time of the inversion 199 procedure, we run the simulation for 30 pacing cycles before measuring the 200 action potential for each iteration in the continuation method used for in-201 version (see Section 2.2.2). However, the parameter changes between each 202 iteration of the continuation method are expected to be quite small and we 203 update the initial conditions between each of the 20 iterations. Therefore, 204 the final solutions of the inversion are expected to be quite close to the steady 205 state solutions for the applied parameters. This has also been confirmed by 206 numerical experiments. For example, the APD90 value for the inversion of 207 data for 10 µM of quinidine was 321.4 ms using 30 pacing cycles from the 208 states saved in the second-to-last continuation iteration. When the simula-209 tion was allowed to continue for 500 pacing cycles, the APD90 value changed 210 by only 0.6% to 323.4 ms. 211

Inversion of action potential measurements 212
The inversion of data of SQT1 hiPSC-CMs exposed to various drugs is per-213 formed using the inversion procedure described in [23]. In the inversion, the 214 adjustment factors λ Kr , λ CaL , λ Na , λ K1 , λ f , λ NaCa , λ NaK , λ bCl , and λ B for 215 the currents I Kr , I CaL , I Na , I K1 , I f , I NaCa , I NaK , I bCl , and the intracellular 216 calcium buffers are treated as free parameters (see (8) and [23]). In addition, 217 we assume that the drugs affect the I Kr , I CaL or I Na currents for quinidine, 218 ajmaline and mexiletine and the I Kr , I f or I Na currents for ivabradine. There-219 fore, we introduce the drug parameters ε Kr , ε f and ε Na for ivabradine and 220 ε Kr , ε CaL , and ε Na for the remaining drugs. In the inversion procedure, we wish to minimize a cost function of the form 223 where D d represent each of the drug doses included in the data set (including 224 the control case, D 0 = 0), H i represent the different cost function terms 225 included in the cost function, w d,i are weights for each of the cost function 226 terms and doses and H reg is a regularization term (see (15) below). The cost 227 function consists of terms of the form where R i (λ, ε, D d ) is a biomarker computed from the solution of the model 229 specified by λ and ε for the drug dose D d , and R * i (D d ) is the corresponding 230 measured biomarker of the data we are trying to invert.

Considered biomarkers
We consider each of the biomarkers, R i , in-232 cluded in the data sets from [9,11], i.e., APD50, APD90, APA, dvdt and 233 RMP. These biomarkers are illustrated in the left panel of Figure 3. Here, 234 the maximal upstroke velocity, dvdt (in mV/ms), is defined as the maximum 235 value of the derivative of the membrane potential with respect to time, the 236 resting membrane potential, RMP (in mV), is defined as the minimum value 237 of the membrane potential, the action potential amplitude, APA (in mV), 238 is defined as the difference between maximum and minimum values of the 239 membrane potential, and the APD50 and APD90 values are defined as the 240 time (in ms) from the time point of the maximum upstroke velocity to the 241 time point where the membrane potential first reaches a value below 50% or 242 90%, respectively, of the action potential amplitude.

243
Cost function weights We use the weight 3 for H APA , the weight 6 for 244 H APD90 and the weight 1 for the remaining terms. In addition, the weights 245 for the control case, w 0,j , and the weight for H APD90 for the largest dose are 246 multiplied by the total number of doses included in the data set. The term 247 H RMP is only included for the drug quinidine because it was not specified in 248 the data set for the remaining drugs [9,11]. 249 Regularization term The regularization term, H reg (λ, ε) is defined as Here, the first term is defined to make the inversion procedure avoid solu-251 tions with unrealistically large perturbations of some of the currents, and the 252 second term is defined to make the inversion select small drug perturbations 253 if small and large perturbations result in almost indistinguishable solutions 254 (D is the median of the considered drug doses of the data set). In our com-255 putations, we set the weights for the terms to w λ = 1 and w ε = 0.001. We 256 let the set S ε consist of all the currents we assume may be affected by the 257 drug and the set S λ consist of the I Kr and I CaL currents. The action potential durations (in ms) at 90% and 50% repolarization (APD90 and APD50, respectively), the maximal upstroke velocity (dvdt, in mV/ms), the action potential amplitude (APA, in mV) and the resting membrane potential (RMP, in mV). Right: Illustration of the QTp interval (in ms) computed from a pseudo-ECG waveform. The QTp interval is defined as the time from the onset of the QRS complex to the peak of the T-wave.

Minimization method 259
We apply the continuation-based minimization procedure from [23] to min-260 imize the cost function (13). We use 20 continuation iterations with 100 or 261 200 randomly chosen initial guesses for the first fifteen and the last five iter-262 ations, respectively. The initial guesses for λ are chosen within 20% above or 263 below the optimal values from the previous iteration, and the initial guesses 264 for ε m are chosen within [ε m−1 /5, 5ε m−1 ], where ε m−1 is the optimal ε from 265 the previous iteration. From these initial guesses we run 30 or 60 iterations of 266 the Nelder-Mead algorithm [42] for the first fifteen and the last five iterations, 267 respectively. In order to represent the WT and SQT1 I Kr currents, parameters of the model for the I Kr open probability are fitted to data of WT and SQT1 I Kr currents from [14]. In this data, it is revealed that the steady state inactivation of the current is shifted towards higher potentials in the SQT1 (N588K) case compared to the WT case (see the right panel of Figure 4). Therefore, we assume that the mutation can be represented in the model by an adjusted model for the steady state value of the inactivation gate, x Kr2 . To find a model matching the data as well as possible, we introduce the six free parameters c 1 − c 6 in the steady state activation and inactivation gates of the I Kr current on the form: See the Supplementary Information for a description of the full I Kr model 270 formulation. To find optimal values of the six parameters c 1 − c 6 , we run 271 simulations of the voltage clamp protocol applied in [14]. That is, we first fix 272 the membrane potential at -80 mV and run a simulation to steady state, then 273 we increase the membrane potential to a value between -50 mV and 100 mV 274 and compare the current after 2 seconds of simulation to the corresponding 275 current reported in [14] (see left panel of Figure 4). In addition, we compare 276 the steady state values of the inactivation gate x Kr2 directly to the steady 277 state inactivation measurements reported in [14].

278
In the simulations trying to adjust the I Kr model to measurements from 279 [14], we adjust the temperature and intracellular and extracellular potassium 280 concentrations in the model to match those reported for the experiments, that 281 is T = 37 • C, [K + ] i = 130 mM, [K + ] e = 4 mM. 282 We find the optimal parameters c = (c 1 , ..., c 6 ) by minimizing a cost 283 function of the form using the Nelder-Mead algorithm [42].  In order to estimate drug effects on adult QT intervals, we apply a simple 290 pseudo-ECG calculation applied in a number of earlier studies (e.g., [43, 44, 291 45, 46, 47]). The calculation follows a two-step procedure.

292
Step 1 The membrane potential and membrane currents along a strand of 293 cylindrical cells are computed using the cable equation as explained i detail 294 in [35]. In this step, cells are connected to each other by gap junctions, each 295 cell is assumed to be isopotential, and the extracellular potential is assumed 296 to be constant in space. The membrane potential v k in each cell k is modeled 297 by where r is the cell radius, L is the cell length, R CG is the ratio between the 299 capacitive and the geometrical cell areas (see e.g., [35]), C m is the specific 300 cell capacitance and σ k−1/2 i is the averaged intracellular conductivity between 301 cells k − 1 and k. This averaged intracellular conductivity is given by where R myo is the myoplasmic resistance and R g is the gap junction resis-303 tance. Furthermore, I k ion is the sum of the ionic current densities of the cell, 304 modeled by the base model (i.e., j I j in (1)).

305
The ODE system defined by (20) is solved using the ode15s solver in 306 Matlab, and the solution is used to compute the membrane currents 307 I k,n m = πr 2 σ originating from each cell i at each time point n.

308
Step 2 The extracellular potential originating from the cell strand is com-309 puted for each time step n using the so-called point-source approximation 310 (see e.g., [43,48]) where σ e is the extracellular conductivity and |r−r k | is the Euclidean distance 312 between the point at which we are measuring the extracellular potential, r, 313 and the center of cell number k, r k .

314
Parameter values The parameters of the pseudo-ECG simulations are 315 based on the parameters of earlier computations of pseudo-ECGs (e.g., [43, 316 46, 35]). We consider a cell stand of 100 cells of length L = 150 µm and 317 radius r = 10 µm. We let the cell strand exhibit transmural heterogeneity 318 of ion channel density, with an endocardial region consisting of the first 25 319 cells, a midmyocardial region for the next 35 cells and an epicardial region 320 in the last 40 cells. The default adult base model defines the ion channel 321 density in the epicardial region. In the endocardial region, the I to and I Ks 322 channel densities are reduced to 1% and 31%, respectively, compared to the 323 epicardial region. Similarly, in the midmyocardial region, the I to and I Ks 324 channel densities are reduced to 85% and 11%, respectively, (compared to 325 the epicardial region).

326
Furthermore, we set C m = 1 µF/cm 2 , R CG = 2, R myo = 0.15 kΩcm, and 327 R g = 0.0015 kΩcm 2 , resulting in a conduction velocity along the cell strand of 328 approximately 52 cm/s, close to experimentally measured conduction veloc-329 ities from literature (∼50 cm/s [49]). Stimulation is applied for the first two 330 cells of the endocardial region, and the extracellular potential is measured 331 2 cm from the end of the epicardial region. Because we in this study only are 332 interested in the time course of the ECG and not the amplitude in mV, we 333 do not assign a specific value for σ e and report the computed pseudo-ECG 334 without numbers on the u e -axis.

335
Definition of the QT interval After computing a pseudo-ECG using the 336 described approach, the QT interval is computed from the ECG waveforms. 337 The QT interval is often defined as the time from the start of the QRS 338 interval to the end of the T-wave (see Figure 3). However, in the study [15], 339 which we will use for comparison to the computational results, an alternative 340 definition of the so-called QTp interval is applied, defined as the time from 341 the start of the QRS complex to the peak of the T-wave, as illustrated in the 342 right panel of Figure 3. We will therefore use this definition in this study. 343 Note that for the SQT1 case, the computed T-wave is inverted to have the 344 opposite sign of the QRS complex (see Figure 5). In that case, we define the 345 peak of the T-wave as the time point when the minimum, i.e. the maximum 346 absolute deviation of the T-wave, is reached. 347 3 Results

348
In this section, we describe the results of our applications of the above men-349 tioned methods. First, we set up the default hiPSC-CM and adult base 350 models for WT and SQT1 based on measurements of WT and SQT1 I Kr 351 currents from [14]. Next, we apply the inversion procedure to identify the 352 effect of four drugs on individual ion channels based on action potential 353 measurements of SQT1 hiPSC-CMs from [9,11]. Finally, we estimate the 354 v (mV)  [14]. Right panel: Comparison of the steady state inactivation gate in the WT and SQT1 models of I Kr and steady state inactivation data from [14]. In the SQT1 case, the inactivation is shifted towards higher values of the membrane potential.
corresponding drug effect for adult patients with short QT syndrome using 355 these identified drug effects.  Figure 4, 363 the fitted model solutions of I Kr in the WT and SQT1 cases are compared 364 to the data from [14]. In the right panel, we observe that the inactivation is 365 shifted towards higher values of the membrane potential in the SQT1 case. 366 The WT and SQT1 versions of the I Kr model fitted to measurements 367 from [14] are inserted into the base model formulation. The parameters, 368 λ hiPSC , for the hiPSC-CM versions of the model are then fitted to action 369 potential measurements of WT and SQT1 hiPSC-CMs from [9]. In addition, The computed APD90 values for hiPSC-CMs are compared to data from [9], and the computed QTp intervals for adults are compared to data from [15].
about the QTp interval for adults with and without the SQT1 mutation 372 from [15]. These adjustments are made by hand-tuning the parameters of 373 the default hiPSC-CM and adult versions of the base model from [23].

374
The action potentials of the resulting models are plotted in the upper 375 panel of Figure 5, along with the computed pseudo-ECG for the adult case. 376 The pseudo-ECG is computed using the approach described in Section 2.4. 377 Note that the only difference between the WT and SQT1 versions of the mod-378 els is a difference in the steady state inactivation gate of the I Kr current (see 379 (18)), and that the remaining parameters of the models are the same in the 380 WT and SQT1 cases. Note also that the adjustment of the I Kr steady state 381 inactivation gate is exactly the same in the hiPSC-CM and adult versions of 382 the model (following assumption 2 in Figure 2), but that the geometry of the 383 cell and the density of different types of membrane and intracellular proteins 384 are different between the models for hiPSC-CMs and adult CMs (following 385 assumption 1 in Figure 2).

386
In the lower panel of Figure 5, we report the APD90 values of the WT and 387 SQT1 models in the hiPSC-CM and adult cases. In addition, we report the 388 QTp interval computed from the pseudo-ECG (see Figure 3). In the hiPSC-389 CM case, we compare the APD90 values of the model to the corresponding 390 values reported in [9]. The APD90 value is reduced from 321 ms in the WT 391 case to 188 ms in the SQT1 case, in good agreement with the values reported 392 in [9]. In the adult case, the APD90 value is similarly reduced from 272 ms 393 to 188 ms. Furthermore, the QTp interval is reduced from 292 ms in the WT 394 case to 200 ms in the SQT1 case, close to reported QTp values from [15]. In [9,11], data are reported from measurements of action potential biomark-399 ers for hiPSC-CMs with the SQT1 mutation N588K exposed to drugs at-400 tempting to make the properties of the SQT1 cells more similar to WT cells. 401 In this paper, we wish to use this data to estimate the corresponding drug 402 responses for adult SQT1 patients. In order to make these predictions, we 403 first need to identify the effect of the drugs on individual ion channels in the 404 hiPSC-CM case, based on the data provided in [9,11]. These drug effects 405 are predicted using the inversion procedure described in Section 2.2 using the 406 SQT1 hiPSC-CM model from Figure 5 as a starting point for the inversion. 407 Figure 6 shows how well the biomarkers of the fitted SQT1 hiPSC-CM 408 model match the corresponding biomarkers reported in the data. We observe 409 that the biomarkers of the fitted models are quite similar to the values of the 410 data. In particular, the model seem to capture the APD90 increase resulting 411 from the drugs quite well. In the lower panel of Figure 7, the action potentials 412 of the SQT1 hiPSC-CM models fitted to the data are plotted for the control 413 case and for each of the drug doses, along with the action potential of the 414 default base model for WT hiPSC-CM. We observe that each of the drugs 415 increase the action potential duration, but not enough to fully recapture the 416 WT hiPSC-CM action potential length.

417
The upper panel of Figure 7 reports the drug effects identified by the 418 inversion procedure in the form of ε values. Recall that ε i is defined as 419 1/IC j 50 , where IC j 50 is the drug concentration that blocks the current j by 420 50%. In other words, a high ε value corresponds to a significant drug effect. 421 We observe that the inversion procedure predicts that all four drugs have a 422  Figure 6: Comparison of measured and computed biomarkers. We consider the biomarkers reported in the data from [9,11] (green) and computed from the fitted SQT1 hiPSC-CM models returned by the inversion procedure described in Section 2.2 (purple) for the drugs quinidine, ivabradine, ajmaline and mexiletine. Note that the definition of each of the biomarkers are illustrated in Figure 3 and that RMP data are only included for the quinidine case. Data are shown as the mean ± SEM (standard error of the mean). In [11], measured drug effects on SQT1 I Kr currents are reported for the 424 maximum dose of each of the considered drugs. In order to get an impression 425 of the accuracy of the IC 50 values estimated by the inversion procedure, we 426 compare the block percentages reported in [11] to the corresponding block 427 percentages estimated by the computational procedure. The results are dis-428 played in Figure 8. We observe that the estimated block percentage is very 429 similar to the measured valued for the drugs quinidine and ivabradine. For 430 ajmaline and mexiletine, the computational procedure seems to predict a 431 slightly larger block percentage than what was found in the measurements. 432 However, all identified block percentages are within the standard error of the 433 mean reported for the measurements. Following assumption 3 in Figure 2, we assume that the effect of a drug on 437 a single ion channel is the same for hiPSC-CMs and adult CMs. Using this 438 assumption, we can estimate the effect of the drugs quinidine, ivabradine, 439 ajmaline and mexiletine for adult SQT1 CMs by inserting the ε values iden-440 tified in Figure 7 into the adult SQT1 base model. In Figure 9,  Figure 10: Comparison of measured and estimated adult QTp intervals. We compare the QTp intervals (see Figure 3) estimated by the computational procedure and the measured QTp intervals from [15] for different heart rates. We consider the WT case, the SQT1 case and the SQT1 case with cells exposed to 6.5 µM quinidine. The drug effect estimated by the procedure is based on measured drug effects for SQT1 hiPSC-CMs from [9]. ventricular action potentials for the SQT1 case with no drug and with each of 443 the considered drug doses present. In addition, we show the action potential 444 in the WT case. The lower panel of Figure 9 displays the adult QTp in-445 tervals computed from pseudo-ECG simulations as explained in Section 2.4. 446 We observe that the drugs increase the QTp interval for the SQT1 case to 447 be more similar to the WT QTp interval. For 10 µM of quinidine, the QTp 448 interval is estimated to increase from 200 ms to 242 ms. Similarly, the QTp 449 interval is predicted to increase to 232 ms for 10 µM of ivabradine, 235 ms 450 for 30 µM of ajmaline and 227 ms for 100 µM of mexiletine. For comparison, 451 the WT QTp interval is 292 ms.

452
The effect of the drug quinidine on the QTp interval of adult SQT1 pa-453 tients is reported for a number of different heart rates in [15]. In Figure 10, 454 we compare the predicted adult drug response of quinidine from our com-455 putational procedure to the measured responses in [15]. In the study [15], 456 mean serum concentration was reported as 2.1 mg/L. Applying the molecular 457 weight of quinidine of 324.4 g/mol [50], we consider the drug concentration 458 6.5 µM in the computational model. Moreover, we use the ε values iden-459 tified based on measurements of SQT1 hiPSC-CMs from [9] (see Figure 7). 460 The QTp intervals are computed for different pacing frequencies using the 461 pseudo-ECG approach described in Section 2.4. In Figure 10, we observe 462 that the QTp intervals predicted by the model for different heart rates fit 463 well with the measured data in both the WT case, in the SQT1 case and in 464 the SQT1 case with quinidine. 465

466
In this paper, we have applied the computational procedure introduced in 467 [22, 23] to data of hiPSC-CMs derived from an SQT1 patient, found in [9,11]. 468 In this procedure, we first identified the drug response on individual ion 469 channels based on measurements of drug effects for SQT1 hiPSC-CMs and 470 then estimated the drug effect for an adult SQT1 patient by inserting these 471 single channel drug effects into a model for adult cells. Short QT syndrome is characterized by an abnormally short duration of the 475 QT interval, and several different genetic mutations have proven to give rise 476 to different subtypes of the syndrome [1,2,3]. SQT1-SQT3 are caused by 477 gain of function mutations in the potassium channels responsible for the I Kr , 478 I Ks and I K1 currents, respectively, while SQT4-SQT6 are associated with 479 loss of function mutations in genes encoding calcium channels. Furthermore, 480 SQT7 and SQT8 have been characterized by a loss of function mutation in a 481 gene encoding sodium channels and a mutation in a gene encoding the cardiac 482 Cl/HCO 3 exchangers AE3, respectively [51]. Because of the different effects 483 of different mutations, it is likely that pharmacological treatment of SQT 484 syndrome must be tailored to the individual subtypes. Furthermore, in order 485 to study the effects of mutations computationally, different mathematical 486 representations are needed for the mechanistic causes of the shortening.

487
In this study, we have considered the SQT1 subtype, caused by a mutation 488 (N588K) leading to increased I Kr currents. Based on measurements of WT 489 and mutated I Kr currents [14], we represented the mutation by shifting the 490 inactivation of the channels towards more positive potentials (see (18) and 491 Figure 4). This computational representation of the mutation follows the 492 same lines as [52,53,46], which used similar adjustments of the model for 493 the I Kr , I Ks and I K1 open probabilities to represent the SQT1, SQT2 and 494 SQT3 subtypes, respectively.

495
In the base model used in this study (see the Supplementary Information), 496 as well as in other commonly used action potential models (e.g., [26, 27, 497 54, 30]), the open probability of voltage-gated ion channels are governed by 498 gating variables following the Hodgkin-Huxley formalism. However, more 499 detailed and versatile models for the open probability have been introduced 500 in the form of more complicated Markov models (see e.g., [55,56] for a 501 comparison of the two formalisms). Such Markov models have been suggested 502 to give a more realistic representation of both the effect of mutations and 503 the effect of drugs (see e.g., [57,58,59]), and was applied in computational 504 studies of SQT1 in [45,60,61] and SQT3 in [62]. On the other hand, complex 505 Markov models are associated with a significant increase in the number of 506 free parameters to determine in the models, and we have therefore chosen to 507 use a more simplified Hodgkin-Huxley formalism in this study.

508
The WT and SQT1 versions of the I Kr current were inserted into the base 509 models for hiPSC-CMs and adult CMs to define WT and SQT1 versions of 510 these action potential models. The resulting computed action potentials 511 displayed a clear action potential shortening for the SQT1 case both for 512 hiPSC-CMs and adult CMs (see Figure 5). The computed pseudo-ECGs for 513 the adult case also displayed a clear reduction in the duration of the QT 514 interval for SQT1. Furthermore, the APD90 values of the hiPSC-CM model 515 and the QT intervals of the adult model appeared to fit well with measured 516 data for WT and SQT1 from [9,15]. Because of the possibly lethal risks associated with SQT syndrome, there is 521 an urgent need for efficient therapies, including pharmacological treatment. 522 In [9,11], a number of drug candidates were tested on hiPSC-CMs generated 523 from an SQT1 patient, and the drugs quinidine, ivabradine, ajmailine and 524 mexiletine seemed to be the most promising for increasing the action poten-525 tial duration of the SQT1 hiPSC-CMs. Other drugs shown to block WT I Kr 526 (sotalol, ranolazine, flecainide, amiodarone) did not seem to work as well for 527 hiPSC-CMs affected by the SQT1 mutation. This highlights the need for 528 mutation-specific investigations of drug effects. Here we have tried to extend 529 such investigations computationally, and have used the measurements from 530 [9,11] to estimate the effect of quinidine, ivabradine, ajmailine and mexile-531 tine for adult SQT1 patients using the computational procedure introduced 532 in [22,23]. In [22,23], drug effects were identified in an inversion procedure based on 535 optical measurements of the membrane potential and cytosolic calcium con-536 centration of hiPSC-CMs in a microphysiological system [63]. In the present 537 study, we have used a different type of data in the inversion procedure, i.e., 538 membrane potential measurements from patch-clamp recordings. This alter-539 native data type could have both advantages and disadvantages compared 540 to the optical measurements used in [22,23]. A disadvantage is that we 541 only have measurements of the membrane potential and not any information 542 about the calcium transient. Information about the calcium transient has 543 been shown to potentially improve the identification of parameters in math-544 ematical action potential models [22,64,65]. However, an advantage of the 545 patch-clamp recordings is that we can obtain information about the actual 546 value of the membrane potential instead of just relative pixel intensities from 547 optical measurements of voltage sensitive dyes. Therefore, we can get infor-548 mation about the resting membrane potential, the action potential amplitude 549 and the upstroke velocity, which could all be useful for the identification of 550 parameters. For example, the upstroke velocity is difficult to obtain from 551 optical measurements of the action potential, which makes it is very hard to 552 identify drug effects on the fast sodium current, I Na , unless additional data, 553 e.g. extracellular measurements are included (see e.g., [66]).

554
With information about the upstroke velocity, it should be possible to 555 identify I Na . However, full identification of all model parameters based on 556 membrane potential measurements is very challenging, and in some cases 557 several combinations of parameter values result in virtually identical action 558 potentials (see e.g., [67,65]). In order to investigate which currents were 559 identifiable in the SQT1 hiPSC-CM base model, we applied the singular 560 value decomposition approach from [68]. The result of this analysis is given 561 in the Supplementary Information. In short, the analysis reveals that the I Na , 562 I CaL and I Kr currents are the most identifiable currents of the model. These 563 are the three main currents investigated for drug effects in our study. In 564 addition, we estimate drug effects on I f for the I f -blocking drug ivabradine. 565 However, we note that the identifiability of this current is not very high and 566 that the estimated drug effect for I f consequently is associated with a large 567 degree of uncertainty. Using the computational inversion procedure described in Section 2.2, we 570 estimated drug effects on major ion channels based on action potential mea-571 surements of SQT1 hiPSC-CMs. For quinidine, ajmaline and mexiletine, we 572 considered drug effects on the currents I Kr , I CaL and I Na , and for ivabradine, 573 we considered drug effects on the currents I Kr , I f and I Na , because these 574 currents have been shown to be affected by ivabradine [69,70,71].

575
Each of the considered drugs had a considerable effect of the I Kr current, 576 with predicted IC 50 values of 8.3 µM for quinidine, 17 µM for ivabradine, 577 44 µM for ajmaline and 223 µM for mexiletine (see Figure 7). The effect 578 of a drug can be very different for WT I Kr channels compared to that of 579 mutated I Kr channels (see e.g., [72,73]). For instance, in [72], the IC 50 value 580 of quinidine was found to be 0.6 µM for WT I Kr channels and 2.2 µM for 581 the SQT1 mutation N588K. Consequently, it is not reasonable to compare 582 the IC 50 values identified for the SQT1 case to IC 50 values found in WT I Kr 583 experiments. However, we observe that the identified IC 50 value for quinidine 584 is quite close to the IC 50 value reported for the SQT1 case in [72] (8.3 µM vs 585 2.2 µM). Furthermore, the predicted block percentage for 10 µM of quinidine 586 is very close to the measured block percentage reported in [11] (see Figure 8). 587 The block percentages for the remaining drugs are also quite close to the 588 percentages reported in [11], indicating reasonable identified IC 50 values for 589 SQT1 I Kr . For ivabradine, the identified IC 50 value for I Kr (17 µM) is also 590 in good agreement with the IC 50 value of 10.3 µM reported for I Kr affected 591 by the SQT1 (N588K) mutation in [69].

592
In addition, the inversion procedure identified a significant block of I Na 593 for the drug mexiletine (see Figure 7). This is reasonable, since mexiletine is 594 known as an I Na blocker [74]. However, the identified IC 50 value of 210 µM is 595 larger than IC 50 values found in literature (7.6 µM -43 µM) [75,76]. More-596 over, ivabradine is an I f blocker with an IC 50 value found to be about 2 µM 597 in [70]. The computational procedure is not able to correctly identify this 598 effect, and underestimates the drug effect to an IC 50 value of about 50 µM. 599 This discrepancy is not unexpected, since the identifiability of I f is predicted 600 to be low (see Figure S1 in the Supplementary Information). The IC 50 value 601 for I Na block of ivabradine, on the other hand, was estimated to 47 µM, in 602 good agreement with the measured IC 50 value of 30 µM reported in [71]. 603 Possible effects of ajmaline and quinidine on I CaL and I Na and of mexiletine 604 on I CaL are underestimated by the inversion procedure compared to the IC 50 605 values reported in [76]. This could indicate that the biomarkers considered 606 in the inversion procedure (see Figures 3 and 6) are not enough to properly 607 characterize these less pronounced drug effects or that the applied base model 608 formulation or simplified drug modeling are not sufficient to represent these 609 effects. It also could be that the regularization term used for the ε values can 610 cause underestimation due to its pressure on finding the minimum effect in 611 the chosen fit. It should also be noted that IC 50 values identified in different 612 experiments tend to vary considerably. Using the single channel drug effects identified from membrane potential 615 measurements of hiPSC-CM, we estimated the drug effects in the adult case 616 based on the three assumptions explained in Figure 2. First, we assumed 617 that the density of different types of ion channels (and other membrane and 618 intracellular proteins) differ between hiPSC-CMs and adult CMs, but that 619 the function of a single channel is the same for adult CMs and hiPSC-CMs. 620 This means that a set of parameters representing the density of different types 621 of proteins, e.g. ion channel conductances (see (4)), are different between 622 the hiPSC-CM and adult CM models, but that the models for, e.g., the 623 open probability of the channels are the same in both versions of the model. 624 Second, we assumed that a mutation affects the function of an ion channel in 625 exactly the same manner for hiPSC-CMs and adult CMs. Therefore, we used 626 the same adjustment of the open probability of the I Kr current (see Figure 4) 627 to represent the SQT1 mutation N588K for both hiPSC-CMs and adult CMs 628 (see Figure 5). Third, we assumed that a drug affects a single ion channel 629 in exactly the same manner for hiPSC-CMs and adult CMs, and we could 630 therefore estimate the drug effect in the adult case based on the drug effect 631 identified for hiPSC-CMs.

632
This modeling framework was used in [22,23] to predict side effects of 633 drugs for adult WT CMs based on measurements of WT hiPSC-CMs. How-634 ever, the procedure can also take the effect of mutations into account (as 635 specified in the second assumption of Figure 2), and in the present study, we 636 have used the procedure to estimate drug effects on the ventricular action 637 potential and the QT interval of adult SQT1 patients based on measurements 638 of hiPSC-CMs generated from an SQT1 patient. Because the effect of one of the drugs considered in this study has been 641 investigated for adult SQT1 patients in the form of measured QT intervals 642 of patients' ECGs, we wished to be able to compare the drug responses 643 predicted by our computational procedure to these measured QT intervals. 644 However, the mathematical base model from [23] only represents the adult 645 ventricular action potential and cannot be directly used to compute the adult 646 QT interval. On the other hand, it is generally recognized that there is a 647 clear link between the duration of the ventricular action potential and the QT 648 interval of a patient's ECG. For example, in [77], a simple linear relationship 649 between the APD90 value and the QT interval was observed. Conversely, 650 in the computational study [78], it was shown that there was not a direct 651 linear realtionship between drug effects on the action potential duration and 652 drug effects on the QT interval. For example, drugs blocking the I Na current 653 gave rise to prolongations of the computed QT intervals, even though no 654 prolongation was observed for the action potential duration. In this study, 655 we have chosen to predict the QT interval by computing a pseudo-ECG using 656 a very simple mathematical representation, as applied in several previous 657 studies of drug effects on the QT interval (e.g., [43,46,47]). The applied 658 method appeared to give rise to reasonable pseudo-ECG waveforms and QT 659 intervals (see e.g., Figure 5), but it is a very simplified representation of the 660 dynamics underlying the ECG, and more realistic spatial modeling, possibly 661 extended to whole-heart simulations including a surrounding torso like in 662 [78,79,80,81], could potentially improve the reliability of the computed QT 663 intervals. The estimated effects of the four drugs quinidine, ivabradine, ajmaline and 666 mexiletine on the ventricular action potential and QT interval of adult SQT1 667 patients revealed that all four drugs are expected to lead to a significant in-668 crease in the action potential duration and the QT interval of SQT1 patients 669 (see Figure 9). It is generally hard to validate the drug responses predicted 670 by the computational procedure because data of drug responses for adult 671 SQT1 patients would be required. However, for the drug quinidine, such 672 experiments are reported in [15] in the form of measured QT intervals for a 673 number of different heart rates. These measurements have been conducted 674 for SQT1 patients without pharmacological treatment, SQT1 patients us-675 ing the drug quinidine and for healthy adults without the SQT1 mutation 676 [15]. We found that the effect of quinidine for adult patients estimated by 677 the computational procedure fitted well with these measured drug responses 678 (see Figure 10), indicating that the applied computational procedure shows 679 promise for reliably predicting adult drug responses based on measurements 680 of hiPSC-CMs.

682
In this paper, we have used a computational procedure to predict the ef-683 fect of the drugs quinidine, ivabradine, ajmaline and mexilitine for clinical 684 biomarkers of adult SQT1 patients, based on measured drug responses for 685 hiPSC-CMs generated from an SQT1 patient. The procedure consists of iden-686 tifying drug effects on individual ion channels based on membrane potential 687 measurements of drug effects for hiPSC-CMs and then mapping this effect 688 up to the adult case based on the assumption of functional invariance of ion 689 channels during maturation. All four drugs were predicted to lead to signif-690 icant increases in the ventricular action potential duration and QT interval. 691 Furthermore, the effect of quinidine predicted by the computational proce-692 dure appeared to be in good agreement with measured effects of the drug. 693 Consequently, we conclude that the computational procedure applied in this 694 study could be a useful tool for estimating adult drug responses based on 695 measurements of hiPSC-CMs, also as new measurements become available, 696 including measurements of alternative drugs or mutations.