#PAGE_PARAMS# #ADS_HEAD_SCRIPTS# #MICRODATA#

Intermolecular interactions play a role in the distribution and transport of charged contrast agents in a cartilage model


Authors: Jenny Algotsson aff001;  Peter Jönsson aff001;  Jan Forsman aff002;  Daniel Topgaard aff001;  Olle Söderman aff001
Authors place of work: Division of Physical Chemistry, Lund University, Lund, Sweden aff001;  Division of Theoretical Chemistry, Lund University, Lund, Sweden aff002
Published in the journal: PLoS ONE 14(10)
Category: Research Article
doi: https://doi.org/10.1371/journal.pone.0215047

Summary

The transport and distribution of charged molecules in polyelectrolyte solutions are of both fundamental and practical importance. A practical example, which is the specific subject addressed in the present paper, is the transport and distribution of charged species into cartilage. The charged species could be a contrast agent or a drug molecule involved in diagnosis or treatment of the widespread degenerative disease osteoarthritis, which leads to degradation of articular cartilage. Associated scientific issues include the rate of transport and the equilibrium concentrations of the charged species in the cartilage and the synovial fluid. To address these questions, we present results from magnetic resonance micro-imaging experiments on a model system of articular cartilage. The experiments yield temporally and spatially resolved data on the transport of a negatively charged contrast agent (charge = -2), used in medical examinations of cartilage, into a polyelectrolyte solution, which is designed to capture the electrostatic interactions in cartilage. Also presented is a theoretical analysis of the transport where the relevant differential equations are solved using finite element techniques as well as treated with approximate analytical expressions. In the analysis, non-ideal effects are included in the treatment of the mobile species in the system. This is made possible by using results from previous Monte Carlo simulations. The results demonstrate the importance of taking non-idealities into account when data from measurements of transport of charged solutes in a system with fixed charges from biological polyelectrolytes are analyzed.

Keywords:

magnetic resonance imaging – cartilage – Mass diffusivity – Finite element analysis – Solute transport – Convection – Electrostatics

Introduction

The function of articular cartilage is to provide low friction and wear resistance, as well as to distribute load in load-bearing joints, for instance, in hips and knees [1]. The optimal functioning of the joint depends on the structure, composition as well as the integrity of the extracellular matrix in cartilage, which is mainly composed of water, collagen and proteoglycans [1, 2]. The load resistance property of cartilage is largely dependent on the proteoglycan content, originating from the glycosaminoglycan (GAG) side chains of the proteoglycan [3]. The GAG is highly negatively charged on account of its carboxyl and sulfate groups that are ionized under physiological conditions. These negative charges are fixed within the extracellular matrix and thus the expression fixed charge density (FCD) of the cartilage is used. Please note that the value of FCD is negative.

Osteoarthritis is a degenerative disease that affects the function of cartilage. The disease is a common cause of disability and it is therefore of importance to understand its progression and to develop efficient treatments. In early stages of the disease, a loss of GAG occurs [4] and one method used to monitor the resulting decrease in the FCD is delayed gadolinium-enhanced magnetic resonance imaging of cartilage (dGEMRIC) [5]. This method utilizes magnetic resonance imaging (MRI) together with an intravenous injection of a contrast agent. The contrast agent commonly used is Gd(DTPA)2− which is negatively charged and contains the paramagnetic ion Gd3+, which has a concentration dependent influence on the spin-lattice relaxation time, T1, of the 1H of the water [6]. The dGEMRIC method is based on the fact that the Gd(DTPA)2− complex will, due to its negative charge, distribute in an inverse relation to the GAG content of the cartilage. This means that where there is a high concentration of GAG, indicating a healthy cartilage, there will be a low concentration of contrast agent, and vice versa. Therefore, a measure of the concentration of the Gd(DTPA)2−, and consequently a proxy for the GAG content of cartilage, is obtained by measuring T1 of the water. In passing, we note that the dGEMRIC method is used both in vitro and in vivo [5, 712]. The method has found most of its applications in in vitro studies and this will most likely continue to be the most important application.

dGEMRIC is dependent on the transport of Gd(DTPA)2− from the site of the injection to the joint of interest, and, subsequently, on the diffusion controlled transport into cartilage, since cartilage lacks direct blood supply. It is therefore important to understand the physiochemical mechanisms that govern the transport and partitioning of Gd(DTPA)2− in cartilage. The transport of the contrast agent will, for example, be affected by the chemical potential difference between the cartilage and the synovial fluid and, in this respect, the FCD of the cartilage plays an important role. Furthermore, the transport will be affected by the amount of collagen and glycosaminoglycan on account of steric effects [1317]. This is an example of the general and important problem of diffusion of small molecules in a complex micro-heterogeneous environment, where crowding and excluded volume effects as well as electrostatic effects are important. The results presented here are therefore applicable to a wider class of related problems than just issues dealing with cartilage.

A previously developed magnetic resonance micro-imaging (μMRI) setup for the investigation of features related to the dGEMRIC method in a model system of cartilage was used to investigate the questions raised above under controlled conditions [18]. The μMRI setup gives spatially and temporally resolved information on the transport of Gd(DTPA)2− in the model system, which is designed with the goal of investigating the effects of the electrostatic interactions in cartilage. To this end, the cartilage is represented by a polyelectrolyte solution and the synovial fluid is represented by a NaCl solution. A semipermeable membrane is used to separate the polyelectrolyte from the salt solution. The experimental results are compared with numerical solutions obtained using the finite element method (FEM) of the relevant equations describing the transport of ionic species in our system. Since previous studies have shown that ions in cartilage do not behave ideally [1921], the influence of intermolecular interactions are included in the theoretical work presented in this paper. This is made possible by comparing the obtained data with previous Monte Carlo simulations on the model system [21]. Furthermore, the low concentration of Gd(DTPA)2− makes it possible to describe the transport with an approximate analytic solution for a system without boundaries.

In summary, our research goal is to address the transport of charged molecules in complex micro-heterogeneous systems and to develop methods where intermolecular interactions are accounted for.

Model system

The model system is outlined in Fig 1A–1C. There is a polyelectrolyte solution at z < 0 and a salt solution at z > 0, separated by a semipermeable membrane at z = 0. The charged polyelectrolyte cannot pass through the membrane into the salt solution, which leads to a difference in electric potential on the two sides of the membrane, which in turn will be compensated for by a concentration difference of ions. The following scenario with respect to time, t, is investigated:

  • At t < 0, there is equilibrium of Na+ and Cl ions between a salt solution and a solution containing a negatively charged polyelectrolyte.

  • At t = 0, a small amount of Gd(DTPA)2− (see Fig 1D for the chemical structure of the complex) is added to the salt solution. The concentration of Gd(DTPA)2− is much lower than the concentration of the other charged species in the system.

  • At t > 0, the Gd(DTPA)2− will diffuse into the polyelectrolyte solution.

Presentation of the experimental and theoretical model systems.
Fig. 1. Presentation of the experimental and theoretical model systems.
(A) The model system at different times. The polymer is negatively charged. At t < 0 there is thermodynamic equilibrium between the two solutions. At time t = 0, Gd(DTPA)2− is added to the salt solution and is allowed to flow over the semipermeable membrane into the polyelectrolyte solution. In the system c Gd ( DTPA ) 2 − , p / s ≪ c Cl − , p / s, c Na + , p / s, FCD (the subscripts p and s denote polyelectrolyte- and salt-solution, respectively). Moreover, because the polymer is negatively charged, the following will be true: c Gd ( DTPA ) 2 − , p < c Gd ( DTPA ) 2 − , s, c Cl − , p < c Cl − , s and c Na + , p > c Na + , s. (B) Schematic drawing of the sample holder used in the μMRI experiments. (C) Cross section of the sample holder in the x-y plane with indicated dimensions. (D) Chemical structure of Gd(DTPA)2− [6]. (E) Chemical structure of carboxymethylcellulose (CMC).

The concentration difference of Gd(DTPA)2− between the solutions at t > 0 will lead to transport of the Gd(DTPA)2− into the polyelectrolyte solution. Initially, the difference is large and leads to fast transport by mutual diffusion from the salt solution into the polyelectrolyte solution (a flow of Na+ and Cl ions will also occur to maintain the electro-neutrality in the system, which is minor compared to the total concentration of Na+ and Cl ions). As equilibrium is approached, the rate of transport decreases. Since the transport is only driven by diffusion, it will take a significant length of time to reach equilibrium. Recall that the root mean square displacement is proportional to the square root of time for diffusive motion.

Theoretical background

The basic equations used in the interpretation of the data are presented in this section. The full derivation of the relations used in the Results and Discussion section is presented in the accompanying S1 and S2 Appendices.

The diffusive flux of an ionic species i, Ji, occuring during step 3 above, can be described by the generalised Fick’s first law (in one dimension):


where Di is the diffusion coefficient and ci the concentration of ion species i, R is the gas constant and T the temperature. Furthermore, μi is the (electro)chemical potential of ion species i, which can, in turn, be expressed by the following equation:

where μ i θ is the standard chemical potential, μ i corr a correction factor for the chemical potential (resulting from intermolecular interactions), and zi the charge of ion species i. F is the Faraday constant and ϕ the electric potential. If the expression for μi in Eq 2 is inserted into Eq 1, we obtain:

Moreover, the change of ion concentration with time, t, can be obtained from the equation of continuity:

The electric potential, ϕ, is related to the charge distribution through Poisson’s equation:


where ϵ0 is the permittivity of vacuum, ϵr is the dielectric constant and FCD is the fixed charge distribution, where FCD = 0 in the salt solution (z > 0) and non-zero in the polyelectrolyte solution (z < 0). In the present system, the concentration of Gd(DTPA)2− is 2 to 3 orders of magnitude lower than the concentration of Na+ and Cl and, as a consequence, Eq 5 is reduced to:

and the electric potential will thus be approximately constant for t > 0 on account of the fact that the relative changes in concentrations of Na+ and Cl are small at these times.

To be able to solve Eqs 3 and 4, μ i corr also needs to be known. It can be estimated by using results from previously made Monte Carlo simulations performed on a system developed to capture the electrostatic interactions in the cartilage/synovial fluid system at equilibrium [21]. In the simulations, the following holds for the two separate solutions, if the electric potential in the salt solution is set as the reference:



where μ i ex is the the excess chemical potential for ion i (related to the activity coefficient, γi, according to: μ i ex = R Tlnγ i) that is resulting from intermolecular interactions (i.e. long-ranged electrostatic interactions and repulsion because of steric overlap of particles), Δϕ is the difference in the electric potential between the two solutions, also known as Donnan potential, and subscripts s and p denote salt solution and polyelectrolyte solution, respectively. The simulations represent a situation of equilibrium, i.e. μi,s = μi,p. The simulations also show that μ i ex could approximately be seen as constant in the system over the steps 1, 2 and 3 above since, the simulations show that at low concentration of Gd(DTPA)2−, the intermolecular interactions between Gd(DTPA)2− and Gd(DTPA)2− ions can be ignored. In the simulations in Ref. [21], the studied concentration of Gd(DTPA)2− was 1 mM.

Material and methods

Materials

Sodium carboxymethylcellulose (CMC) (average molecular weight = 90 kDa, 99.5% purity, Sigma Aldrich; molecular structure is given in Fig 1E) was used without further purification. The degree of substitution of the used batch of the polyelectrolyte has previously been determined to 0.9 [18]. Magnevist® (0.5 mmol gadolinium/mL, Bayer Healthcare) was used as the Gd(DTPA)2− source. The molecular structure of Gd(DTPA)2− is given in Fig 1D. NaCl (99.5% purity, Merck) and NaOH (analytical grade) were used as received. All solutions were prepared with MilliQ purified H2O. Regenerated cellulose acetate membranes with a cutoff of 5 kDa were obtained from Millipore. This cutoff was previously shown not to affect the transfer rate of Gd(DTPA)2− significantly and at the same time only allow for a negligible amount of CMC to pass through the membrane [18]. For dialysis of the polyelectrolyte solutions, Spectra/Por® Biotech Cellulose Ester (CE) Dialysis membrane (dialysis cutoff of 100-500 Da) was used.

Sample preparation

The pH of the salt solution (150 mM NaCl) was adjusted to between 8 and 9 by means of drop-wise addition of 1 M NaOH. This pH value is slightly higher than the generally quoted value of synovial fluid of around 7.5. In this pH range, effectively all of the carboxylic groups on CMC are charged.

Magnevist® was diluted by the salt solution to a Gd(DTPA)2− concentration of 0.36 mM. Before use, the pH was adjusted to the same value as in the salt solution.

When preparing the polyelectrolyte solutions, the CMC powder was dissolved in the salt solution (150 mM NaCl). Before use, the polyelectrolyte solutions (7 mL) were dialyzed for 24 h against 1.5 L 150 mM NaCl solution. A control experiment showed that this time was enough to reach ionic equilibrium between the polyelectrolyte solution and the salt solution (for more information, see S3 Appendix). The salt solution (the dialyzate) was gently stirred during the dialysis and the pH of the dialyzate was kept between 8 and 9 by periodically making drop-wise addition of 1 M NaOH. The carbon content in the polyelectrolyte solutions was determined by total organic carbon (TOC) analysis after the dialysis. The FCD was then calculated using as input the degree of substitution of the polyelectrolyte. To obtain further information with regard to the ionic content, the sodium content was determined in the solutions by inductively coupled plasma atomic emission spectroscopy (ICP-AES) analysis. Subsequently, the Cl concentration could be determined by applying electroneutrality conditions in the solutions. The densities of the solutions were determined by a precision density meter (DMA5000, Anton-Paar). The polyelectrolyte concentration of the polyelectrolyte solutions was estimated to 1.7, 2.2 and 2.5 wt%, respectively.

A small piece of the regenerated cellulose acetate membrane was punched out to fit into the sample holder (see Fig 1B and 1C). Before use, the membrane was, according to instructions from the manufacturer, put in 150 mM NaCl solution (with the glossy side up) to wash out glycerol and to wet the pores. The solution was changed three times during 1 h and thereafter, the membrane was left in the solution for at least 24 h before use.

Experiments

The bottom section of the sample holder (constructed from polyether ether ketone (PEEK) plastic and a 10 mm NMR tube; see Fig 1B and 1C) was filled with ∼160 μL CMC solution (added by weight). Thereafter, the membrane was carefully placed on top of the polyelectrolyte solution with the glossy side down. Subsequently, the top part of the sample holder (including the NMR tube) was assembled with the bottom part and 205 μL of the 150 mM NaCl solution was added on top of the membrane.

2.1 mL of the 0.36 mM Gd(DTPA)2 solution was added to the top compartment about 30 minutes after the addition of the salt solution (see above). The first T1 measurement (for description of T1 measurements, see below) was started within a few minutes after the addition of Gd(DTPA)2−. Due to rapid initial changes of the Gd(DTPA)2− concentration, the initial T1 measurement (for each sample) gave unphysical concentration profiles in the salt solution and these data points were not used in the subsequent analysis. The sample holders were kept in the spectrometer for the initial 22 h and T1 measurements were performed continuously. Subsequently, the sample holders were removed from the spectrometer. For acquisition of subsequent data, the sample holder was re-inserted into the spectrometer for each measurement. Data at 3 different positions were collected in each measurement (one set for the whole polyelectrolyte solution, one set simultaneously showing the polyelectrolyte and salt solutions, and one set for the salt solution close to the membrane). This measurement procedure was carried out once a day. Between the measurements, the sample holders were stored at 298 K, protected from light. The samples were not stirred during the experiments and, in order to minimize agitation, the sample holders were carefully handled between measurements. Three different concentrations of CMC were studied and for each concentration, two independent experiments were performed.

μMRI measurements

Spatially resolved T1 measurements were recorded at 298 K on an 11.7 T Bruker Avance II 500 spectrometer equipped with a Bruker MIC-5 microimaging probe having a maximum gradient strength of 3 T m-1 and a 10 mm saddle-coil radio-frequency insert. The measurements were performed using a pulse sequence containing an inversion recovery block followed by a rapid acquisition with relaxation enhancement (RARE) block [22] for the read out of the images. The recovery time delay, τ, was incremented logarithmically in 32 steps from 0.0017 s to 15 s, while the images were taken up with an 11.1 × 4.8 mm field-of-view (z × y), 128 × 48 acquisition matrix size and a 0.2 mm slice thickness (x). This resulted in a measurement time of about 10 minutes per image.

Data processing

The images were reconstructed to a 128 × 48 matrix size, which gives a spatial resolution of 87 × 100 μm, and subjected to 0.2 mm Gaussian smoothing. For each voxel, the T1 values were obtained by fitting the image intensity, I, to following equation:


where I0 is the maximum image intensity, τ is the recovery time delay, and A is a constant, which is less than 2 if the pulse repetition time is less than around 5 times T1 or the inversion pulse is imperfect. Errors representing one standard deviation were obtained by performing Monte Carlo error estimations as described in [23].

The obtained T1 values were then converted to the concentration of Gd(DTPA)2−, c G . 2 −, using the expression:


where R G . 2 − is the relaxivity of Gd(DTPA)2−, T 1 G . 2 − and T1pre are the values of T1 with and without Gd(DTPA)2−, respectively. R G . 2 − was determined in all used solutions (the salt solution and the three different polyelectrolyte solutions with different concentrations of CMC) by measuring T1 for different concentrations of Gd(DTPA)2−. The R G . 2 − value in the salt solution was 3.84 ± 0.03 (mM s)-1. The obtained values for the polyelectrolyte solutions are collected in Table 1. Errors in R G . 2 − are estimated from the uncertainties in the obtained T1 values.

Tab. 1. Parameter values.
Parameter values.
Determined values for the different sets of experiments.

The three different images that were taken at the same time point (see above) were visually merged to obtain a larger field-of-view of the system.

Assuming that there is no convection in the polyelectrolyte solution, the diffusion coefficients of Gd(DTPA)2− in the polyelectrolyte solutions were obtained from fitting the relevant equations to the concentration profiles in the polyelectrolyte solution (see S1 Appendix for explanation of the procedure). The data points used in the fit were obtained in a time interval where the influence of the membrane is low and the Gd(DTPA)2− has not yet reached the bottom of the polyelectrolyte compartment, i.e. 2.5-6 h.

All data processing was performed with in-house written Matlab® code.

Numerical methods

Numerical calculations were carried out using COMSOL Multiphysics® 4.4, a program that solves partial differential equations using the finite element method (FEM). The FEM simulations were carried out for a one-dimensional geometry with the same dimensions as in Fig 1C, except that the salt solution was truncated at z = dsalt to account for thermally driven convection in the salt solution. At this distance, the concentration in the salt solution in the FEM simulations was set equal to the experimental value at the same position. The distance dsalt was varied between 1 to 5 mm for the different experiments, depending on whether the contribution from convection was high or low, respectively. This was determined by comparing the experimental concentration profiles in the salt solution with those obtained from the FEM simulations without convection and choosing dsalt at a value where the experimental and the calculated curves in the salt solution are equal. Moreover, we assume homogenous properties in each solution, such that values of the diffusion coefficients are constant in each solution and the chemical potential correction factors are, except close to the interface, also constant in each solution.

Data from all experiments (two sets of data for each polyelectrolyte concentration) were used in the simulations. First, the time-independent Eq 4 and Poisson’s equation (Eq 6) were solved for the Na+ and Cl ions as well as for the electric potential. In order to avoid a singularity in the change in chemical potential at z = 0, we introduced a smooth change in μ i corr around z = 0 according to:


where Δ μ i corr = μ i , p corr - μ i , s corr (in what follows, the second indices p, and s denote polyelectrolyte and salt solutions, respectively). Choosing a different value than 1 nm for the range over which the chemical potential change has negligible effect on the macroscopic difference in electric potential between the polyelectrolyte and salt solution, Δϕ, and thus according to Equation B in S1 Appendix does not affect the obtained value for the macroscopic step in concentration of Gd(DTPA)2− at z = 0. The obtained values of c Na +, c Cl − and ϕ were subsequently used to calculate the distribution of Gd(DTPA)2− as a function of time by solving Eqs 3 and 4 for Gd(DTPA)2−. The FEM simulations were run for each experiment between 80 h and 260 h.

The boundary conditions used are summarized in Table 2. The initial concentrations were c Na + , s = 150mM, c Cl − , s = 150mM, c G . 2 − , s = 0.33mM and c G . 2 − , p = 0mM. The temperature was 298.15 K, ϵr = 78, and the diffusion coefficients for Na+ and Cl were assigned values of 1.33⋅10-9 m2/s and 2.03⋅10-9 m2/s, respectively [24]. It should be noted that the obtained results are independent of the chosen diffusion coefficients, since we are assuming a situation of steady state for the Na+ and Cl ions. The value for the diffusion coefficient for Gd(DTPA)2− used in the salt solution was obtained from a study of molecules similar to Gd(DTPA)2− and was taken to be 4.5⋅10-10 m2/s [25]. The diffusion coefficient for Gd(DTPA)2− in the membrane was set to 1.5⋅10-10 m2/s in all FEM simulations, chosen to fit with the experimental data at early measurement times. As noted above, the diffusion coefficients in the polyelectrolyte solution were obtained from analysis of the experimental data (see also below and in the S1 Appendix). Values for Δ μ i corr were determined from the equilibrium concentrations together with equations describing equilibrium concentrations in the framework of a Donnan equilibrium. However, values for Δ μ i corr cannot be determined independently for Na+ and Cl as there are two independent relations and three unknowns (Δ μ Na + corr, Δ μ Cl − corr and ΔμG.2−corr, see Equations D-F in S2 Appendix). To arrive at individual values of Δ μ i corr needed in the analysis of Gd(DTPA)2−, we start by setting the Δ μ Na + corr to values of Δ μ Na + ex obtained from Monte Carlo simulations [21]. The Δ μ Cl − corr values are subsequently obtained from the equilibrium concentrations of Na+ and Cl (see Table 1) and Equations D and E in S2 Appendix. In principle, values for Δ μ G . 2 − corr can then be determined from Equation F in S2 Appendix. However, the equilibrium concentrations of Gd(DTPA)2− are less accurately known, on account of the long time required to reach equilibrium. To improve the accuracy in the obtained value, we have made an initial estimate from Equation F in S2 Appendix and then allowed the equilibrium value of Gd(DTPA)2− to vary in the FEM simulations so that the simulated concentration profiles best fit the experimentally determined ones. The used values are collected in Table 1.

Tab. 2. Boundary conditions.
Boundary conditions.
Boundary conditions used in the finite element method (FEM) simulations.

To investigate the validity of the approximate equations that are derived in the Results and Discussion section, FEM simulations were also carried out with same dimensions as in Fig 1C (not truncated at a distance dsalt as described above) and values taken from Monte Carlo simulations (Ref. [21] and unpublished results, see S3 Fig).

Results and discussion

The outline of this section is as follows. We start by presenting and discussing the experimental data obtained, and then proceed to numerically solve the relevant equations describing the data using FEM. Subsequently, we develop approximate solutions that can be used to illustrate how different parameters affect the distribution of ions in the system. Using examples, we demonstrate under which conditions the approximate relations are valid. Finally, we use the approximate relations to determine the diffusion coefficients of the contrast agent in the polyelectrolyte solution. In the interest of keeping the paper short, we delegate derivations to the S1 and S2 Appendices, to which we refer when needed.

Experimental data

The main experimental results of this study are concentration profiles for Gd(DTPA)2− as a function of time for three different values of FCD. An example is given in Fig 2, which shows the experimentally determined concentrations of Gd(DTPA)2− in a system with the polyelectrolyte solution having a FCD of -108 mM, at different times after addition of Gd(DTPA)2− to the salt solution (results for the duplicate experiment with FCD = -108 mM and the other values of FCD are shown in the S1 Fig). After the addition, Gd(DTPA)2− diffuses from the salt solution to the polyelectrolyte solution. This is a slow process as can be observed in Fig 2. The times needed for the average concentration in the polyelectrolyte solution to reach 50% and 90% of the equilibrium concentration in the 10 mm long polyelectrolyte solution are 40 h and 230 h, respectively. However, there is, as expected, a rapid development of a discontinuity in the concentration at z = 0, although it appears to be gradual in Fig 2 due to the finite resolution of the μMRI method and artifacts due to the presence of the membrane (e.g. uncertain values of relaxivities and T1pre). The FEM simulations make it possible to study the discontinuity closer (see below).

Experimental concentration profiles.
Fig. 2. Experimental concentration profiles.
Concentration profiles of Gd(DTPA)2− at different times in the salt solution and the polyelectrolyte solution with FCD = -108 mM (z < 0). Gd(DTPA)2− was injected in the salt solution at t = 0 and the data was acquired using μMRI. The (small) decrease of the concentration of Gd(DTPA)2− at longer times on the right hand side is due to the finite volume of the salt solution.

In Fig 3A, we show the concentration of Gd(DTPA)2− in the polyelectrolyte solution at a position 4 mm from the membrane as a function of time for a FCD of -108 mM. The concentration of Gd(DTPA)2− will initially increase rapidly in the polyelectrolyte solution, but the magnitude of the increase will gradually decrease with time. The symbols present two sets of data recorded in repeated experiments with the same FCD for the polyelectrolyte solution. The difference in the curves is attributed to a different degree of temperature gradient driven convection in the salt solution at the beginning of the experiments. The contribution from convection can be estimated by analyzing the concentration profiles in the salt solution (see Fig 2A and S1 Fig). Furthermore, convection evens out concentration gradients in the salt solution and therefore, the time to reach equilibrium becomes shorter. On account of the geometry (small radius and height) and the position of the sample holder in the spectrometer (where temperature gradients are small), and the fact that the polyelectrolyte is rather viscous, the transport in the polyelectrolyte solution is well described by pure diffusion.

Compiled results from experiment and FEM simulations.
Fig. 3. Compiled results from experiment and FEM simulations.
(A)The concentration of Gd(DTPA)2− in the polyelectrolyte solution at a position 4 mm from the membrane for two samples with FCD = -108 mM: one with small effects due to convection (circles) and one with noticeable convection at the beginning of the experiment (triangles), as a function of time. The solid lines are the corresponding values from FEM simulations with fitted values of Δ μ G . 2 − corr and the dashed line represents the corresponding results from a simulation with Δ μ i corr = 0 for all ions. (B) Experimental (dashed) and simulated (solid) values for the concentration of Gd(DTPA)2− at different positions along the container, z, for a polyelectrolyte solution with FCD = -108 mM. (C) Simulated values for the concentration of Gd(DTPA)2− at different positions along the container, z, for a polyelectrolyte solution with FCD = -108 mM. The dashed line shows the result from a FEM simulation with fitted value of Δ μ G . 2 − corr and the solid line represents the corresponding results from a simulation with Δ μ i corr = 0 for all ions. (D) The step in concentration at z = 0, Δ c G . 2 − = C p ( 0 − ) G . 2 − − C s ( 0 + ) G . 2 − (see also Fig 4), at different times for three polyelectrolyte solutions with different FCD. The dashed lines correspond to the individual experiments and the solid lines represent the corresponding mean value of the duplicate experiments. The data are taken from FEM simulations.

In the following sections, we will use relevant theoretical expressions in combination with FEM simulations to investigate and rationalize the behavior of the Gd(DTPA)2− with regard to its distribution and transport. The analysis yields physical properties relevant for the system. We note that the approach can be used to predict and understand the transport of ions in general in related systems.

Finite element method simulations

By solving Eqs 3, 4 and 6 numerically, using FEM techniques (as described above in Material and methods), we are able to predict the experimental results in Fig 2. Using as input the data in Table 1 (and values for the diffusion coefficients for Na+ and Cl given above) and the boundary conditions in Table 2, the results are presented as solid lines in Fig 3A and 3B (for additional comparison between experiments and FEM simulations, see S2 Fig). As can be seen, there is a good agreement with the experimental results at all times, and the effects due to convection in the salt solutions are well accounted for. We stress that all relevant parameters have been determined in independent experiments, except for the values of the electric potential and the chemical potential correction factor for the mobile ions, which are in good agreement with the results from independent Monte Carlo simulations (Ref. [21] and unpublished results, see S3 Fig).

One of the key motivations for carrying out this study was to investigate to what extent effects due to non-ideal conditions are important in the transport and equilibrium concentrations of Gd(DTPA)2− in cartilage and models of cartilage. This is one important issue in the dGEMRIC method. To this end, we present in Fig 3A and 3C results of FEM simulations where the difference in the chemical potential correction factors between the two solutions of all species are set to zero. It is clear that the predicted equilibrium concentrations deviate considerably, in agreement with earlier work [5, 18, 21].

The FEM simulations enable us to study the discontinuity in the concentration at z = 0. In Fig 3D, the step in the Gd(DTPA)2− concentration, Δ c Gd ( DTPA ) 2 − (see also Fig 4), is plotted against time. There is a rapid development of a discontinuity in the concentration at z = 0 and the magnitude of the step increases with time. The double experiments for each value of FCD, show that the magnitude of the step size is reproducible. The deviations are due to different degrees of convection in the salt solution. While the steps in concentration vary with time, the ratios between the concentrations of Gd(DTPA)2− in the polyelectrolyte solution and salt solution at z = 0 are independent of time, with values of 0.82, 0.78 and 0.74 at FCD values of -73 mM, -92 mM and -108 mM, respectively.

Concentration and potential profiles from FEM simulations.
Fig. 4. Concentration and potential profiles from FEM simulations.
(A)A concentration profile for Gd(DTPA)2− 5h after addition of Gd(DTPA)2− in the salt solution. (B) The corresponding potential profile. The profiles are taken from FEM simulations with a polyelectrolyte solution with a FCD of -108 mM. The inserts to the right show the profiles around z = 0.

Approximate relations for concentration profiles

Motivated by the need to illustrate the influence of the relevant parameters and to carry out rapid and straightforward calculations of the diffusion and concentration profiles, we here develop approximate expressions describing the transport and concentration profiles of Gd(DTPA)2− in our model system. Our point of departure is to consider an interval 0z ≤ 0+ (see Fig 4), the size of which is such that we can assume that quasi steady-state conditions apply. This is justified by the fact that the gradient in the concentration at z = 0 adjusts much more rapidly than the variation in concentration over macroscopic length scales outside 0z ≤ 0+. In the current experiments, this interval is around 10 nm which is orders of magnitude smaller than the dimensions of the macroscopic container (see Fig 4). Solving Eq 3 with Ji = 0 gives (the reader is referred to S1 Appendix for details):


where C s ( 0 + ) i is the macroscopic concentration of ion i in the salt solution at z = 0+. To obtain approximate expressions for the step sizes in potential and concentration between the salt solution and the polyelectrolyte solution, Eq 12 is linearized (see Equation D in S1 Appendix). Inserted into Eq 6, this gives the following, approximate expression for the electric potential drop between the polyelectrolyte solution and the salt solution:

where Δϕ is the electric potential difference between the polyelectrolyte solution and the salt solutions, CNaCl,s is the concentration of Na+ and Cl in the salt solution (the small difference in concentration due to the presence of Gd(DTPA)2− ions is neglected) and Δ μ i corr is defined above. We note that Ohshima and coworkers have analyzed a system that is similar to ours and present results which are in agreement to those presented here when excess chemical potentials are neglected [26, 27]. By combining Eqs 12 and 13 and linearizing the exponent, we obtain the following approximate expressions for the step in concentration, Δci, between the polyelectrolyte solution and the salt solution:



From Eqs 14 and 15, it is clear that the difference in the chemical potential correction factor between the solutions changes the step in the concentration at the interface with the same amount for both Na+ and Cl, and therefore, as noted above, the individual values of Δ μ i corr for the different ions cannot be determined from measurements of the concentration steps alone. The expressions above are only valid for small changes in the concentration (on account of the linearization), but they are useful in order to illustrate how the different parameters affect the distribution of ions in the system and qualitatively understand the behaviour of Fig 3D. Thus, the concentration of Gd(DTPA)2− in the polyelectrolyte solution will decrease due to the negative value of the FCD of the polyelectrolyte solution. Furthermore, the difference in μ i corr between the two solutions can either increase the concentration step, as for Na+, or decrease the concentration step, as for Cl and Gd(DTPA)2− (cf. Eqs 1416 and Table 1). The step size will increase as Gd(DTPA)2− reaches the bottom wall of the container, resulting in an increase in C s ( 0 + ) G ⋅2 −, which is the reason for the increasing magnitude of the concentration step size observed in Fig 3D. Additionally, Eqs 1416 illustrates the fact that the assumption of Δ μ i corr = 0 for all mobile ions, mentioned above, has a substantial impact on the step size at all times (cf. Fig 3C and Table 1).

We next carry out an approximate analysis of the transport of Gd(DTPA)2−. To this end, we solve Eq 17 on both sides of the salt/polyelectrolyte solution interface at z = 0


which is obtained by combining Eqs 3 and 4 and is valid outside the interval 0 < z < 0+. The term containing the electric potential in Eq 3 can be neglected since the concentration of Na+ and Cl is much larger than the concentration of Gd(DTPA)2− resulting in the electro-kinetic transport of the Gd(DTPA)2− ions being small compared to the diffusional transport. The concentration difference between the salt solution and the polyelectrolyte solution is approximately given by Eq 16, which together with the condition that the flux of ions on both sides of z = 0 should be equal, couples the solutions for the polyelectrolyte and salt solutions. The solution to these equations for an infinite system without convection and membrane is given by (see S1 Appendix for details):


where C s , 0 G . 2 − is the concentration of Gd(DTPA)2− at t = 0 in the salt solution (z < 0) and D G . 2 − , p and D G . 2 − , s are the diffusion coefficients of Gd(DTPA)2− in the polyelectrolyte solution and the salt solution, respectively. The concentrations Cp(0−)G.2− and C s ( 0 + ) G . 2 − are related via Eq 16, where C p ( 0 − ) G ⋅ 2 − = C s ( 0 + ) G . 2 − + Δ c G . 2 −. From the condition that the flux is equal at z = 0 and at z = 0+, the following condition also applies:

where the ratio C p ( 0 − ) G . 2 − / C s ( 0 + ) G . 2 − can be obtained from Eq 16. The expressions in Eqs 18 and 19 are also approximately valid for short times in a bounded system for times on the order of t max < 0.1 L p 2 / D G . 2 − , p where Lp is the length of the polyelectrolyte solution. With Lp = 1 cm and D G . 2 − , p = 4 ⋅ 10 − 10 m2/s, this gives tmax < 7 h, which should be considered as an order of magnitude estimate above which the effect of the boundaries in the system starts to influence the concentration profile significantly. However, the expressions in Eqs 18 and 19 give an estimate of the time scales needed for the Gd(DTPA)2− ions to diffuse in the system. The time it takes for the polyelectrolyte solution to reach a certain average value would then be expected to scale as L p 2.

Approximate analytic expressions vs finite element simulations

To illustrate under which conditions the approximate relations above are valid, we compare them with FEM solutions for the exact (within the model) equations. Fig 5 shows the concentration profiles using the analytical expressions in Eqs 18 and 19 and Equation W in S1 Appendix and how it compares to the simulated values at different times for a low FCD = -50 mM (Fig 5A), for the highest FCD used in this study, i.e. -108 mM (see Fig 5B), and for a higher FCD = -150 mM (see Fig 5C). It can be seen that the approximate analytic expressions reproduce the concentration profiles for values of FCD of -50 mM and -108 mM, respectively, but less so for FCD = -150 mM on account of the linearization of Equation B in S1 Appendix. We note that for short times (less than 2 h), it is non-trivial to include the effect of the membrane on the distribution of ions in the different compartments (see S1 Appendix for further information). The concentration profile in the polyelectrolyte solution starts to deviate considerably from the analytical expressions when the Gd(DTPA)2− ions have reached the wall of the container containing the polyelectrolyte solution (not shown).

Comparison of concentration profiles from FEM simulations and approximate analytic expressions.
Fig. 5. Comparison of concentration profiles from FEM simulations and approximate analytic expressions.
The concentration of Gd(DTPA)2− vs z from FEM simulations (solid lines) and from Eqs 18 and 19 and Equation W in S1 Appendix (dashed lines) for a polyelectrolyte solution with an FCD of (A) -50 mM, (B) -108 mM and (C) -150 mM. Parameters used in the calculations A; D G . 2 − , p = 4 . 2 ⋅ 10 − 10 m2/s, Δ μ Na + corr = − 0 . 0736R T, Δ μ Cl − corr = − 0 . 041R T, Δ μ G . 2 − corr = − 0 . 125R T. B; values in Table 1. C; D G . 2 − , p = 3 . 5 ⋅ 10 − 10 m2/s, Δ μ Na + corr = − 0 . 140R T, Δ μ Cl − corr = − 0 . 095R T, Δ μ G . 2 − corr = − 0 . 318R T.

Diffusion coefficient of contrast agent in the polyelectrolyte solution

From Eq 18, it is possible to obtain values of the diffusion coefficient of Gd(DTPA)2− in the polyelectrolyte solution from the experimentally determined concentration profiles. As outlined in S1 Appendix, one has to take into consideration the fact that Eq 18 applies to a system with infinite boundaries and that there is convection in the salt solution. When this is accounted for, Eq 18 is modified to Equation Z in S1 Appendix. A fit of Equation Z in S1 Appendix to one dataset is given in S1 Appendix, and the obtained diffusion coefficients are presented in Table 1. We note that the values follow the expected dependence on the magnitude of FCD. Thus, the values are 91, 87 and 82% of the value in the salt solution for FCD values of -73, -92 and -108 mM, respectively. As noted in the Introduction, steric effects may affect the diffusion coefficient of the contrast agent, which raises the question whether such specific effects differ between our model system and real cartilage. It is a general observation that electrostatic effects dominate diffusion of charged species in charged colloidal systems of different kinds. Moreover, results in Refs. [1317] indicate that steric effects are relatively minor in real cartilage. Taken together these observations does not prove, but support the notion that the diffusion data for the contrast agents in our model system has bearing to the situation in real cartilage. Finally, we note that the problem of diffusion of charged species through charged colloidal systems is a general problem in a multitude of applications, and we feel that our suggested approach of obtaining such data is one of the main results of this work.

Conclusions

We have combined magnetic resonance imaging results with simulations on a previously developed model system with the aim of studying issues related to early diagnosis and in vitro and in vivo clinical studies of osteoarthritis. The flow of a charged contrast agent from a salt solution into a polyelectrolyte solution has experimentally been determined with spatial and temporal resolutions of 90 μm and 600 seconds, respectively. The relevant equations describing the flow have been solved using Finite Element Methods. The results clearly show that effects due to non-ideality must be taken into account when analyzing flow and equilibrium distributions of divalent species, such as contrast agents, in complex polyelectrolyte systems, such as cartilage. Approximate expressions for the flow have been derived and compared with numerical solutions of full equations in order to assess when the approximations are reasonable. Finally, the derived approximate expressions make it possible to determine the diffusion coefficient for the contrast agent in a polyelectrolyte solution from the experimental data.

Supporting information

S1 Appendix [pdf]
Diffusion-driven transport of ions in an infinite, one-dimensional Donnan system.

S2 Appendix [pdf]
Donnan equilibrium.

S3 Appendix [pdf]
Determination of dialysis time.

S1 Fig [pdf]
Additional experimental concentration profiles.

S2 Fig [pdf]
Additional concentration profiles from experiments and FEM simulations.

S3 Fig [pdf]
Obtained values of from Monte Carlo simulations.


Zdroje

1. Mow VC, Ratcliffe A, Poole AR. Cartilage and diarthrodial joints as paradigms for hierarchical materials and structures. Biomaterials 1992; 13: 67–97. doi: 10.1016/0142-9612(92)90001-5 1550898

2. Flik KR, Verma N, Cole BJ, Bach BR. In Cartilage Repair Strategies; Williams RJ. Ed.; Humana Press: Totowa, NJ, 2007; pp 1–12.

3. Buschmann MD, Grodzinsky AJ. A Molecular Model of Proteoglycan-Associated Electrostatic Forces in Cartilage Mechanics. Journal of Biomechanical Engineering 1995; 117: 179–192. doi: 10.1115/1.2796000 7666655

4. Buckwalter JA, Mankin HJ, Grodzinsky AJ. Articular Cartilage and Osteoarthritis. Instr. Course Lect. 2005; 54: 465–480. 15952258

5. Bashir A, Gray ML, Burstein D. Gd-DTPA2− as a measure of cartilage degradation. Magn. Reson. Med. 1996; 36: 665–673. doi: 10.1002/mrm.1910360504 8916016

6. Caravan P, Ellison JJ, McMurry TJ, Lauffer RB. Gadolinium(III) Chelates as MRI Contrast Agents: Structure, Dynamics, and Applications. Chem. Rev. 1999; 99: 2293–2352. doi: 10.1021/cr980440x 11749483

7. Bashir A, Gray ML, Hartke J, Burstein D. Nondestructive imaging of human cartilage glycosaminoglycan concentration by MRI. Magn. Reson. Med. 1999; 41: 857–865. doi: 10.1002/(sici)1522-2594(199905)41:5<857::aid-mrm1>3.0.co;2-e 10332865

8. Trattnig S, Mlynárik V, Breitenseher M, Huber M, Zembesch A, Rand T, et al. MRI visualization of proteoglycan depletion in articular cartilage via intravenous administration of Gd-DTPA. Magn. Reson. Imaging 1999; 17: 577–583. doi: 10.1016/s0730-725x(98)00215-x 10231184

9. Nieminen MT, Töyräs J, Lassanen MS, Silvennoinen J, Helminen HJ, Jurvelin JS. Prediction of biomechanics proper ties of articular cartilage with quantitive magnetic resonance imaging. Journal of Biomechanics 2004; 37: 321–328. doi: 10.1016/s0021-9290(03)00291-4 14757451

10. Wang N, Chopin E, Xia Y. The effects of mechanical loading and gadolinium concentration on the change of T1 and quantification of glycosaminoglycans in articular cartilage by microscopic MRI. Phys. Med. Biol. 2013; 58: 4535–4547. doi: 10.1088/0031-9155/58/13/4535 23760174

11. Sigurdsson U, Siversson C, Lammentausta E, Svensson J, Tiderius CJ, Dahlberg LE. In vivo transport of Gd-DTPA2− into human meniscus and cartilage assessed with delayed gadolinium-enhanced MRI of cartilage (dGEMRIC). BMC Musculoskeletal Disorders 2014; 15: 226. doi: 10.1186/1471-2474-15-226 25005036

12. Söderman O, Algotsson J, Dahlberg LE, Svensson J. Biophysics and Biochemistry of Cartilage by NMR and MRI; The Royal Society of Chemistry, 2017; pp 176–190.

13. Salo EN, Nissi MJ, Kulmala KAM, Tiitu V, Töyräs J, Nieminen MT. Diffusion of Gd-DTPA2− into articular cartilage. Osteoarthritis and Cartilage 2012; 20: 117–126. doi: 10.1016/j.joca.2011.11.016 22179030

14. Bansal PN, Stewart RC, Entezari V, Snyder BD, Grinstaff MW. Contrast agent electrostatic attraction rather than repulsion to glycosaminoglycan affords a greater uptake ratio and improved quantitative CT imaging in cartilage. Osteoarthritis and Cartilage 2011; 19: 970–976. doi: 10.1016/j.joca.2011.04.004 21549206

15. Kulmala KAM, Karjalainen HM, Kokkonen HT, Tiitu V, Kovanen V, Lammi MJ, et al. Diffusion of ionic and non-ionic contrast agents in articular cartilage with increased cross-linking—Contribution of steric and electrostatic effects. Medical Engineering & Physics 2013; 35: 1415–1420.

16. Arbabi V, Pouran B, Weinans H, Zadpoor AA. Multiphasic modeling of charged solute transport across articular cartilage: Application of multi-zone finite-bath model. Journal of Biomechanics 2016; 49: 1510–1517. doi: 10.1016/j.jbiomech.2016.03.024 27033729

17. Pouran B, Arbabi V, Zadpoor AA, Weinans H. Isolated effects of external bath osmolality, solute concentration, and electric charge on solute transport across articular cartilage. Journal of Biomechanics 2016; 49: 1510–1517.

18. Algotsson J, Forsman J, Topgaard D, Söderman O. Electrostatic interactions are important for the distribution of Gd(DTPA)2− in articular cartilage. Magn. Reson. Med. 2016; 76: 500–509. doi: 10.1002/mrm.25889 26332213

19. Stell G, Joslin CG. The donnan equilibrium: a theoretical study of the effects of interionic forces. Biophys. J. 1986; 50: 855–859. doi: 10.1016/S0006-3495(86)83526-3 19431690

20. Dai H, Potter K, McFarland EW. Determination of ion activity coefficients and fixed charge density in cartilage with 23Na magnetic resonance microscopy. J. Chem. Eng. Data 1996; 41: 970–976. doi: 10.1021/je9600257

21. Algotsson J, Åkesson T, Forsman J. Monte Carlo Simulations of Donnan Equilibrium. Magn. Reson. Med. 2012; 68: 1298–1302. doi: 10.1002/mrm.24409 22890897

22. Hennig J, Nauerth A, Friedburg H. RARE imaging: A fast imaging method for clinical MR. Magn. Reson. Med. 1986; 3: 823–833. doi: 10.1002/mrm.1910030602 3821461

23. Alper JS, Gelb RI. Standard Errors and Confidence Intervals in Nonlinear Regression: Comparison of Monte Carlo and Parametric Statistics. J. Phys. Chem. 1990; 94: 4747–4751. doi: 10.1021/j100374a068

24. Cussler EL. Diffusion: Mass Transfer in Fluid Systems; Cambridge University Press, 2nd edn., 1997.

25. Vander Elst L, Sessoye A, Laurent S, Muller RN. Can the Theoretical Fitting of the Proton-Nuclear-Magnetic-Relaxation-Dispersion (Proton NMRD) Curves of Paramagnetic Complexes Be Improved by Independent Measurements of Their Self-Diffusion Coefficients? Helvetica Chimica Acta. 2005; 88: 574–587. doi: 10.1002/hlca.200590040

26. Ohshima H, Miyajima T. Evaluation of the Donnan Model for Polyelectrolytes Using the Composite Poisson-Boltzmann Equations. Colloid Polym. Sci. 1994; 272: 803–811. doi: 10.1007/BF00652421

27. Ohshima H. The Donnan Potential-Surface Potential Relationship for a Cylindrical Soft Particle in an Electrolyte Solution. J. Colloid Interface Sci. 2008; 323: 313–316. doi: 10.1016/j.jcis.2008.04.027 18502443


Článok vyšiel v časopise

PLOS One


2019 Číslo 10
Najčítanejšie tento týždeň
Najčítanejšie v tomto čísle
Kurzy

Zvýšte si kvalifikáciu online z pohodlia domova

Získaná hemofilie - Povědomí o nemoci a její diagnostika
nový kurz

Eozinofilní granulomatóza s polyangiitidou
Autori: doc. MUDr. Martina Doubková, Ph.D.

Všetky kurzy
Prihlásenie
Zabudnuté heslo

Zadajte e-mailovú adresu, s ktorou ste vytvárali účet. Budú Vám na ňu zasielané informácie k nastaveniu nového hesla.

Prihlásenie

Nemáte účet?  Registrujte sa

#ADS_BOTTOM_SCRIPTS#