# Driven translocation of a semiflexible polymer through a conical channel in
the presence of attractive surface interactions

Andri Sharma<EMAIL_ADDRESS>Rajeev Kapri<EMAIL_ADDRESS>Abhishek Chaudhuri<EMAIL_ADDRESS>Department of Physical
Sciences, Indian Institute of Science Education and Research (IISER) Mohali,
Sector 81, S. A. S. Nagar 140306 Punjab, India

###### Abstract

We study the translocation of a semiflexible polymer through a conical channel
with attractive surface interactions and a driving force which varies
spatially inside the channel. Using the results of the translocation dynamics
of a flexible polymer through an extended channel as control, we first show
that the asymmetric shape of the channel gives rise to non-monotonic features
in the total translocation time as a function of the apex angle of the
channel. The waiting time distributions of individual monomer beads inside the
channel show unique features strongly dependent on the driving force and the
surface interactions. Polymer stiffness results in longer translocation times
for all angles of the channel. Further, non-monotonic features in the
translocation time as a function of the channel angle changes substantially as
the polymer becomes stiffer, which is reflected in the changing features of
the waiting time distributions. We construct a free energy description of the
system incorporating entropic and energetic contributions in the low force
regime to explain the simulation results.

##  I Introduction

The translocation of macromolecules through nanopores and channels has been a
topic of immense interest since the pioneering work by Kasianowicz et al. in
which DNA was successfully translocated through a biological nanopore
Muthukumar (2016). Study of confined polymers and translocation of polymers
through geometrical pores paved a way to study sequencing of DNA and RNA,
which has a lot of applications in the field of medicine and biotechnology.
During the last few decades, a number of experimental and theoretical work has
been performed to understand the dynamics of the transportation of such
biomolecules Palyulin  _et al._ (2014); Branton  _et al._ (2010); Wanunu
(2012); Howorka and Siwy (2009); Keyser (2011); Milchev (2011); Bezrukov  _et
al._ (1994); Kasianowicz  _et al._ (1996); Sung and Park (1996); Park and Sung
(1998); Muthukumar (1999, 2003); Sakaue (2007, 2010); Saito and Sakaue (2012);
Sarabadani  _et al._ (2014); Rowghanian and Grosberg (2011); Ikonen  _et al._
(2012a, 2013, b); Bhattacharya and Binder (2010); Bhattacharya (2013); de Haan
and Slater (2010); Gauthier and Slater (2008a, b); Huang  _et al._ (2014); Luo
_et al._ (2007, 2008); Cohen  _et al._ (2011, 2012a); Sarabadani and Ala-
Nissila (2018); Suhonen and Linna (2018); Katkar and Muthukumar (2018); Kumar
Sharma  _et al._ (2019); Suma and Micheletti (2017); Suma  _et al._ (2020);
Kumar Sharma  _et al._ (2019); Plesa  _et al._ (2016).

The transport of DNA through a nanopore/nanochannel is facilitated by the
application of an electric field across the channel, giving rise to an
electrostatic force which pulls the DNA into the channel. Due to the reduction
in conformational entropy of a polymer in confinement, the geometry of the
nanochannel has been found to strongly affect the dynamics of the
translocation process. For example, in DNA sequencing using biological
nanopores like $\alpha$-hemolysin and MspA, it has been observed that the
shape of channel plays a vital role in determining the ionic current through
it Schneider and Dekker (2012); Sun and Luo (2014); Sun  _et al._ (2018); Wong
and Muthukumar (2010); Meller  _et al._ (2001); Meller (2003); Cohen  _et al._
(2012b); Kumar  _et al._ (2018). Compared to the $\alpha$-hemolysin pore,
where the cylindrical shape of the beta barrel dilutes the ion current
specific to individual nucleotide of the DNA inside the channel, thus making
sequencing difficult, the narrow constriction of the cone shaped MspA allows
better resolution of current signatures corresponding to individual nucleotide
Derrington  _et al._ (2010); Branton  _et al._ (2010); Derrington  _et al._
(2010); Venkatesan and Bashir (2011).

Figure 1: Snapshots from the simulation for various stages of the
translocation process of a semiflexible polymer with $N=64$ beads, through a
conical channel of length $L_{p}=16\sigma$, half apex angle $\alpha$ and
channel entrance width $d$. (a) Equilibrium polymer conformation on the cis
side of the channel with one end fixed at the pore entrance. Subsequent
polymer conformations (b) midway through the translocation process, (c)
showing a possible hairpin formation near the channel exit and (d) when the
polymer has successfully translocated to the trans side of the channel.

In a bid to control the transport dynamics of the DNA inside the nanochannel
to achieve better sequencing, experimental approaches have shifted to the
building of bio-inspired nanopores and nanochannels with similar transport
properties Dekker (2010); Wanunu and Meller (2007); Gershow and Golovchenko
(2007); Luan  _et al._ (2012); Keyser  _et al._ (2006). These artificial
solid-state nanochannels not only allow the control of the shape of the
channel and tunability of its surface properties, but are also stable with
respect to changes in pH, temperature and mechanical oscillations.
Asymmetrical conically shaped nanopores have been found to have excellent
sensing applications, with the nanaopore tip acting as a sensing zone Mathé
_et al._ (2005); Harrell  _et al._ (2006); Lan  _et al._ (2011); Thacker  _et
al._ (2012); Zhou  _et al._ (2017); Chen  _et al._ (2021). Conical
nanocapillaries show ionic current rectification at low salt concentrations.
Further, chemical modifications of the inner surfaces of a conical nanochannel
can lead to reversal in the ion current rectification direction. A combination
of a conical nanofunnel with a cylindrical nanochannel, was shown to
significantly reduce the threshold external voltage required to trigger the
translocation of DNA through the channel Zhou  _et al._ (2017).

Theoretical studies of polymer translocation through conical nanochannels have
shown that the escape of a confined flexible polymer from a conical channel is
a pore driven process and can proceed without an external force Nikoofard  _et
al._ (2013); Nikoofard and Fazli (2015). The asymmetric shape of the channel
ensures that the polymer has a larger entropic penalty near the constriction
which gives rise to an entropic force leading to the escape of the polymer
from the larger opening. Further, the passage time is a non-monotonic function
of the apex angle of the cone for a given length of the channel. Langevin
dynamics studies of flexible polymer translocation through conical channels
have shown that translocation is dependent on channel structure and
interactions of the polymer with the channel Sun  _et al._ (2018); Tu  _et
al._ (2018); Sun  _et al._ (2021). Most biopolymers and proteins are however
semiflexible with an energy cost associated in bending, characterized by the
bending rigidity $\kappa$ of the polymer. The natural question to ask is how
the stiffness of the polymer affects the translocation dynamics. It is with
this motivation that we characterize the translocation dynamics of a
semiflexible polymer through an interactive conical nanochannel.

In this paper, we present a detailed study of translocation dynamics of both
flexible and semiflexible polymers through cone shaped nanochannels in the
presence of a spatially varying external driving force and surface
interactions, using coarse-grained Langevin dynamics simulations. The external
force acting along the channel length mimics the voltage driven translocation
of polymers, higher at narrower regions of the channel and lower near the
larger channel openings. We first report results of the waiting time
distributions for a translocating flexible polymer for varying strengths of
the external force for a flat channel . The total translocation time $\tau$
decreases with increasing force strength as expected. A break up of the total
translocation time into a filling, transfer and escape time provides valuable
insight on the translocation dynamics. As the apex angle $\alpha$ (see Fig.
1(a)) is varied, the translocation time shows non-monotonic features,
consistent with earlier reports Nikoofard  _et al._ (2013); Nikoofard and
Fazli (2015); Sun  _et al._ (2018); Tu  _et al._ (2018); Sun  _et al._ (2021).
We comment on the detailed behavior of the waiting time distribution with
varying force strengths and $\alpha$. We next present a detailed study of the
translocation dynamics of semiflexible polymers with increasing stiffness
through this conical nanochannel. The non-monotonic features observed in
$\tau$ with varying $\alpha$ differs significantly from that of the flexible
polymer. The translocation time is also strongly dependent on the stiffness of
the polymer. Further, with increasing forces, the non-monotonic features
reduce significantly and total translocation time increases with increasing
$\alpha$ for all values of polymer stiffness. We present phase plots of $\tau$
in the $\kappa-\alpha$ plane and mean waiting time distributions to
characterize the translocation dynamics in detail. Finally, we provide free
energy arguments using a quasi-equilibrium approximation valid at low forces
to explain our observations.

The paper is organized as follows. In section 2, we introduce the simulation
model and methods where the governing equations are explained. This section
also sets the notations used in the paper. In section 3, we present the
results and discussion of our simulations for (a) flexible and (b)
semiflexible polymer and provide possible explanation for the observed
behaviour at low forces via a free energy description. In section 4, we
conclude by discussing the importance of our results and possible future
directions to extend the domain of physical relevance of this work.

Figure 2: External force $f_{ext}$ as a function of the length of the channel
along its axis $x$, plotted for different apex angles $\alpha$ at a fixed
$f_{0}=1.0$ according to Eq. 5. With increasing $\alpha$, external force drops
rapidly with $x$ as expected.

##  II Model

To model the conical channel in two dimensions, we consider the walls of the
channel to be made up of two rows of fixed monomer beads of size $\sigma$ and
of length $L_{p}=N^{\prime}\sigma$, where $N^{\prime}$ is the number of
monomers making up a wall of the channel. The walls are symmetric about the
$x-$axis with the apex angle of the channel $\alpha$ defined with respect to
the $x$-axis as shown in the Fig. 1. The diameter of the channel at the cis
side is fixed at $d$ to allow only single monomer entry. To separate the cis
and trans regions of the channel, vertical walls made up of fixed monomers
beads are placed along the $y-$ direction.

The polymer is modeled via a coarse grained bead-spring chain, where non-
bonded monomers interact via a short ranged repulsive Lennard-Jones (rLJ)
potential given by

$U_{\textrm{bead}}(r)=\begin{cases}4\epsilon\Bigg{[}\left(\frac{\sigma}{r}\right)^{12}-\left(\frac{\sigma}{r}\right)^{6}\Bigg{]}+\epsilon,&r<r_{c}\\\
0,&r\geq r_{c}\end{cases}$ (1)

where $\epsilon$ gives the strength of the potential. The above truncated and
shifted LJ potential has cut-off at $r_{c}=2^{\frac{1}{6}}\sigma$. The
interaction between consecutive monomers of the chain is harmonic with the
interaction given by

$U_{\textrm{bond}}=\frac{1}{2}K(r-r_{0})^{2},$ (2)

where $K$ is the spring constant and $r_{0}$ is the equilibrium separation
between consecutive monomers.

Figure 3: (a) Waiting time distribution $w(s)$ of a flexible polymer through a
flat channel ($\alpha=0^{\circ}$) as a function of scaled translocation
coordinate $s/N$, for weak forces $f_{0}=0.1,0.2$. (b) Waiting time
distribution for strong forces : $f_{0}=1.0,2.0$. (Inset) $w(s)$ for force
values, $f_{0}=4.0$.  Figure 4: Variation of different time components
($\tau_{s}$) with external force $f_{0}$ for a flat channel
($\alpha=0^{\circ}$). Here, $\tau=\tau_{f}+\tau_{p}+\tau_{e}$. Figure 5:
Variation of total translocation time $\tau$ with external force $f_{0}$ for a
flat channel ($\alpha=0^{\circ}$). $\tau\sim 1/f_{0}^{\gamma}$ with
$\gamma\approx 1.78$. Figure 6: Total translocation time ($\tau$) of flexible
polymer as a function of apex angle ($\alpha$) for different external forces
$f_{0}$. (Inset) Magnified plot of $\tau$ versus $\alpha$ at larger force
values showing the disappearance of maxima-minima peaks. Figure 7: Passage
time ($\tau_{p}$) and escape time ($\tau_{e}$) contributions to the total
translocation time for a flexible polymer at a low force value, $f_{0}=0.2$.

To model a semiflexible polymer, an additional bending potential is introduced
between consecutive bonds as follows :

$U_{\textrm{bend}}(\theta_{i})=\kappa(1+\cos{\theta_{i}})$ (3)

where $\theta_{i}$ is the angle between the $i^{th}$ and $(i-1)^{th}$ bond
vectors. $\kappa$ quantifies the bending rigidity of the polymer.

The interaction of the polymer beads with the vertical walls are modelled by
the same repulsive Lennard-Jones (rLJ) introduced before for the polymer. The
surface interaction of the conical channel with the polymer is however
attractive, with the attraction between the beads of the polymer and the
channel beads given by the standard Lennard-Jones interaction as :

$U_{\textrm{channel}}(r)=\begin{cases}4\epsilon\Bigg{[}\left(\frac{\sigma}{r}\right)^{12}-\left(\frac{\sigma}{r}\right)^{6}\Bigg{]},&r\leq
r_{c}^{lj}.\\\ 0,&\text{otherwise}\end{cases}$ (4)

where $r_{c}^{lj}=2.5\sigma$.

Figure 8:  Waiting time distributions ($w(s)$) of a flexible polymer for
different $\alpha$ for $f_{0}=0.2$.

In addition to the forces on the polymer due to its interactions with the
channel, it experiences an additional external force in the positive $x-$
direction, modeled as a potential gradient, Fig. 2 :

$f_{ext}=\frac{f_{0}d}{(d+2x\tan\alpha)}$ (5)

where $f_{0}$ is the constant force expected in a linear channel with
$\alpha=0^{\circ}$. Evidently, for a conical channel the force decreases along
its length.

The equation of motion for the position of the $i$th monomer of the polymer is
given by

$m\ddot{\vec{r_{i}}}=-\zeta\dot{\vec{r_{i}}}-\nabla
U_{i}+\vec{f}_{ext}+\vec{\eta}_{i},$ (6)

here m is the mass of the monomer, $U_{i}$ is the net potential faced by the
monomer, $\zeta$ is the friction coefficient and $\vec{\eta}_{i}$ is a
Gaussian random force with
$\langle\eta_{i}(t)\eta_{j}(t^{{}^{\prime}})\rangle=2k_{B}T\zeta\delta_{ij}\delta(t-t^{{}^{\prime}})$
where $k_{B}$ is the Boltzmann’s constant and $T$ is the temperature. The
equations are simulated in LAMMPS using the Verlet update scheme. The scales
of length, energy and mass are set by $\sigma$, $\epsilon$ and $m$
respectively. This sets the time scale as
$\tau_{0}=\sqrt{\frac{m\sigma^{2}}{\epsilon}}$. In these units, we choose
$k_{B}T=1.0$, $\zeta=1.0$, $r_{0}=1.12$, $d=1.25$ and
$K=10^{3}k_{B}T/\sigma^{2}$.$f_{0}$ is varied in the range
$0.1-2\,\,k_{B}T/\sigma$ and $\kappa$ varies in the range $0-8$. The number of
polymer beads $N=64$ and pore beads $N^{\prime}=16$ are fixed in our
simulations. The time step is set as $\Delta t=0.001\tau_{0}$ and all results
presented are averaged over $1500-2000$ independent samples.

We start our simulation by holding the first bead of the polymer at the narrow
opening of the channel while the rest of the polymer segment is allowed to
relax to its equilibrium conformation outside the channel. The translocation
process is tracked by the translocation coordinate $s$. The counter for $s$
starts once a bead of the polymer reaches the trans side after crossing the
channel. The polymer is said to be successfully translocated if $s$ is equal
to the total number of beads $N$ of the polymer and the simulation terminates.
The speed with which the polymer traverses the channel gives a detailed idea
of what is going on inside the channel. To trace the above effect, we focus on
the waiting time $w(s)$, defined as the average time spent by a monomer inside
the channel.

Figure 9: Total translocation time versus $\alpha$ for a semiflexible polymer
for four different bending rigidties $\kappa=1,2,4$ and $8$. In (a) the
variation is shown for $f_{0}=0.2$. (Inset) Magnified plot at lower $\alpha$
for $\kappa=1.0$ and $2.0$. As external force is increased : (b) $f_{0}=0.5$
(c) $f_{0}=1.0$ and (d) $f_{0}=2.0$, $\tau$ increases almost linearly with
$\alpha$.

##  III Results and Discussion

###  III.1 Flexible Polymer

Figure 10: Translocation time ($\tau$) as a function of apex angle ($\alpha$)
separated into passage time ($\tau_{p}$) and escape time ($\tau_{e}$)
components for two different values of rigidity of the polymer (a) $\kappa=2$
and (b) $\kappa=4$ for $f_{0}=0.2$.

We first look at the translocation dynamics of a flexible polymer through a
flat channel $\alpha=0^{\circ}$ which serves as a control. As the strength of
the external force $f_{0}$ acting inside the channel is increased
($f_{ext}=f_{0}$ for a flat channel), the total translocation time ($\tau$) of
the polymer decreases. To explain this behavior, we look at the waiting time
distribution $w(s)$ of the polymer at different values of $f_{0}$ (Fig. 3).
Tension propagation theory Sarabadani  _et al._ (2014); Ikonen  _et al._
(2012a, b); Bhattacharya (2013) accurately predicts the waiting time
distribution for a slit (channel of unit size). $w(s)$ is expected to rise
with $s$ as the tension due to the external force reaches the end of the
polymer in the cis side. After that the system is in the tail retraction stage
and $w(s)$ starts decreasing with $s$. In our system, where the pore is
extended and is attractive in nature, we observe several interesting features
in $w(s)$ which are different from that of the slit. These arise due to the
interplay between surface interactions and the external drive.

For the initial part of the translocation when the polymer enters the pore,
the attractive nature of the pore and the external force, sucks the polymer
inside the pore. For smaller forces $(f_{0}=0.1,0.2)$ (Fig. 3(a)), $w(s)$
rises with $s$ as the tension slowly reaches towards the end. However, as
$f_{0}$ increases (see Fig. 3(b) with $f_{0}=1.0,2.0$), the combined effect of
the attractive pore and the external force pulls the initial monomers quicker
inside the pore, thus reducing their waiting times. For the lower forces (Fig.
3(a)), we observe a flat regime where the $w(s)$ does not change with
increasing $s$. Subsequently, once the tension has propagated to the end, the
tail retraction part sets in and $w(s)$ starts to fall. Dramatically however,
$w(s)$ shows a sharp rise for the end monomers of the polymer. This feature is
true for all external force values. As the number of monomers left in the
channel becomes lesser than the channel length and keeps decreasing, the net
external force on the polymer decreases. Therefore, the attractive nature of
the channel becomes a dominant factor for the end monomers resulting in
increased waiting times.

The behavior of the waiting times as explained above can be better understood
if we divide the total translocation time as $\tau=\tau_{f}+\tau_{p}+\tau_{e}$
where (i) $\tau_{f}$ is the initial filling time, the time taken by the first
monomer of the polymer to reach the exit without returning back into the
channel; (ii) $\tau_{p}$ is the transfer time, the time taken from the exit of
the first monomer into the trans side to the entry of the last monomer from
the cis side; and (iii) $\tau_{e}$ is the escape time, the time between the
entry of the last monomer in the channel and its escape to the trans-side Luo
_et al._ (2008); Kumar  _et al._ (2018). $\tau_{f}$ depends strongly on the
external drive and the surface interactions which sucks the polymer inside.
$\tau_{p}$ depends on surface interactions and entropy of the polymer segment
outside the channel on the trans side and less on the external force since
these monomers see a constant force inside the channel. $\tau_{e}$ is strongly
influenced by the attractive interactions which hold back the polymer as it
exits the channel. In Fig. 4, we have plotted the three translocation times as
a function of the force. For low forces, $\tau_{f}$ is significantly smaller
than $\tau_{p}$ and $\tau_{e}$, suggesting that translocation is dominated by
the transfer and escape dynamics. This explains the initial rise in $w(s)$ for
small $s$ as successive monomers spent longer times inside the channel. At
larger forces, $\tau_{f}$ becomes comparable to the transfer and escape times.
As explained earlier, the combined effect of surface interactions and external
force, pulls the monomer rapidly inside the pore. Therefore, the initial
monomers show a dip in the waiting times.

In Fig. 5, we have plotted the total translocation time of the polymer as a
function of the external force. The translocation time falls faster than
$1/f_{0}$ for larger forces. This behavior is reminiscent of the effect of
attractive pore-polymer interactions which provides an additional pull on the
monomer beads as they enter the pore. Such behavior of $\tau$ with $f_{0}$ for
driven polymer translocation through extended channels have been reported
earlier Luo  _et al._ (2007).

We next look at the situation for the conical channel ($\alpha\neq 0$). Since
more polymer conformations are possible in the extended part of the channel
compared to the constricted region, this shape asymmetry leads to a force of
entropic origin which can facilitate movement of the polymer from the cis to
the trans end. The shape asymmetry therefore influences the translocation
process significantly in addition to the effects of the attractive surface
interactions of the channel. Further the external force is space dependent
(Fig. 2) and decreases as we move from the cis to the trans side of the
channel. This leads to a global decrease in the total translocation time
($\tau$) when compared to a flat channel. Additionally, $\tau$ displays
distinct non-monotonic features as $\alpha$ is varied (Fig. 6). We first
analyze this non-monotonic behavior for very low forces (say $f_{0}=0.1$).

As $\alpha$ starts to increase, $\tau$ decreases as translocation is
facilitated by the additional entropic drive when compared to
$\alpha=0^{\circ}$ which drives the polymer inside the conical pore faster
compared to a flat channel. This continues till a certain $\alpha=\alpha_{1}$
(see Fig. 6; for $f_{0}=0.1$, $\alpha_{1}\approx 1.5^{\circ}$). Beyond this
apex angle, $\tau$ starts increasing with $\alpha$. Note that with increasing
$\alpha$, the attractive interactions along the length of the channel weaken.
This lowering of the pull on the polymer near the entrance leads to longer
translocation times. It is again useful to split the total translocation time
as was done for the flat channel. However, unlike a flat channel, a conical
channel has a larger exit which results in folding of the polymer back into
the channel (see Fig. 1(c)). Therefore, for such channels it becomes
convenient to split the translocation time into a passage time and an escape
time, $\tau=\tau_{p}+\tau_{e}$. Passage time $\tau_{p}$ for a conical channel
is defined as the time taken from the time of entrance of the first bead of
the polymer to the time the last bead of the polymer reaches the entrance of
the pore from the cis side. The escape time $\tau_{e}$ is defined similar to
that of the flat channel as the time taken after $\tau_{p}$ until all the
polymer beads have exited the pore from the trans side. As we see, from Fig.
7, the non-monotonicity in $\tau$ for smaller $\alpha$ is completely governed
by $\tau_{p}$. The weakening of the surface interactions near the cis side is
what dominates the behavior.

This increase in $\tau$ with $\alpha$ continues till it reaches
$\alpha=\alpha_{2}$. Beyond $\alpha_{2}$, $\tau$ starts to decrease with
$\alpha$. This is a combined effect of the entropic drive and surface
interactions. Due to a larger pore exit, a pore which is completely filled
with the polymer, experiences a stronger entropic drive. Further, the escape
time ($\tau_{e}$) also decreases as the attractive interactions near the exit
becomes negligible and cannot hold back the polymer. With increasing $\alpha$,
we see a curious behaviour. At high $\alpha$, $\tau_{e}$ starts to increase
and becomes larger than $\tau_{p}$ at some $\alpha$. This crossover can be
explained by noting that the exit is now large enough to allow folding of the
polymer leading to hairpin formations. Escape of the polymer becomes the rate
limiting step. Note that beyond a value of $\alpha$, $\tau_{p}$ essentially
saturates as the forces that drive the entry and subsequent passage towards
the exit do not change.

At higher forces ($f_{0}>0.4$; see inset of Fig. 6), $\tau$ increase with
$\alpha$. To understand this behavior, note that for higher forces, the
surface interactions do not play a major role. The force along the channel
axis is larger for smaller $\alpha$. This ensures that at large force values
the translocation process is extremely fast. As $\alpha$ increases, the force
along the channel axis decreases (see Fig. 2) and translocation becomes
slower.

We further verify our findings for the translocation time by looking at the
waiting time distribution of the monomer beads (Fig. 8). At a force value of
$f_{0}=0.2$, we note that as $\alpha$ is increased from $1.5^{\circ}$ to
$3.5^{\circ}$, $w(s)$ for the initial beads increases dramatically. This is a
result of weakening attractive interactions near the cis side as the pore
broadens. Note that for the end beads, $w(s)$ for $\alpha=1.5^{\circ}$ rises
sharply because of the attractive interactions at the trans side which prevent
exit of the polymer. This effect is reduced sharply for $\alpha=3.5^{\circ}$.
As $\alpha$ is increased further, $w(s)$ reduces due to weakening attractive
interactions near the exit and a greater entropic drive. At even higher
$\alpha=10^{\circ}$, we start seeing the effects of the folding in of the
polymer leading to increased $w(s)$ for all monomer beads when compared to
$w(s)$ at lower $\alpha$ values. The dip in $w(s)$ for the last few monomers
at large values of $\alpha$ is easily explained. With surface interactions
lowered significantly and large part of the polymer already outside the pore
leading to a higher entropy, the end monomers are sucked out extremely fast.

Figure 11: Phase plot of total translocation time $\tau$ in the
$\kappa-\alpha$ plane for two forces (a) $f_{0}$=0.2. and (b) $f_{0}=1.0$.
Figure 12: Waiting time distributions for various $\alpha$ at (a) $\kappa=2$,
and (b) $\kappa=4$ for $f_{0}=0.2$. The distributions show distinct
differences as the bending rigidity of the polymer is increased.

###  III.2 Semiflexible Polymer

We now discuss the translocation dynamics of a driven semiflexible polymer
through the conical channel. We first look at the total translocation time as
function of the apex angle for different values of polymer rigidity and
external driving force. In Fig. 9(a), we have plotted $\tau$ as function of
$\alpha$ for a low external force value ($f_{0}=0.2$). The variation of the
translocation time with apex angle still shows a non-monotonic feature as
observed in the case of a flexible polymer. However, unlike that of the
flexible polymer, rigid polymers do not show a sharp secondary peak in $\tau$.
The total translocation time shows an initial decrease with increasing
$\alpha$ indicating that as the apex angle increases, it becomes favorable for
the polymer to exit the channel. This is a result of the enhanced entropic
drive that we discussed for the flexible case. Beyond a critical value of
$\alpha$, $\tau$ starts to increase again. Note that the critical $\alpha$
shifts to larger values of apex angle and to higher values of $\tau$ as the
rigidity of the polymer is increased (see Fig. 9(a)). This indicates that
translocation becomes more difficult with increasing stiffness of the polymer.
This feature is consistent with earlier results of driven translocation of
semiflexible polymers through extended channels Adhikari and Bhattacharya
(2013); Kumar  _et al._ (2018); Sarabadani  _et al._ (2017). The largely
monotonic increase in $\tau$ with increasing $\alpha$ is different from that
of the flexible polymer.

The effect of rigidity can be explained by looking at the multiple stages of
the translocation process. At lower values of rigidity (say $\kappa=2$, see
Fig. 10(a)), $\tau_{p}$ plays a more significant role for most $\alpha$
values. This shows that the total translocation time is dictated by the
filling of the pore as in the flexible case. However, with increasing rigidity
(say $\kappa=4$, see Fig. 10(b)), $\tau_{e}$ becomes the rate limiting
process. The increased influence of the escape time dynamics happens because
the rigid polymer which at this stage fills the channel, keeps deflecting
between the two walls of the channel. This is a behavior typical of
semiflexible polymers with increased rigidity Odijk (1983). The attractive
interactions of the two walls keeps the polymer pinned inside the channel for
long times, increasing the escape time $\tau_{e}$. The monotonic increase of
$\tau$ with $\alpha$ as opposed to the appearance of a second peak for
flexible polymers in the conical channel, can also be explained from the
dominance of the escape time dynamics. The competition of the attractive
interactions with the entropic drive which results in the lowering of
translocation time beyond $\alpha_{2}$ for flexible polymers is no longer the
rate limiting step. Rather, escape time dynamics which are dominated by the
attractive interactions of the channel walls with the escaping beads
influences the translocation process.

As the external force is increased, the non-monotonicity observed in the
variation of $\tau$ with $\alpha$ starts decreasing (Fig. 9(b-d)). At very
high forces, $f_{0}\geq 1.0$ (see Fig. 9(d)), the total translocation time
increases monotonically with the apex angle for all rigidities of the polymer.
In Fig. 11 we construct phase plots for $\tau$ in the $\kappa-\alpha$ plane
for various values of $f_{0}$. For small forces $f_{0}=0.2$, $\tau$ is maximum
for stiffer chains and smaller apex angles (see Fig. 11(a),
$\kappa>3,\,\,0^{\circ}\leq\alpha\leq 1^{\circ}$). For larger forces
$f_{0}=1.0$, the phase plot shows a largely uniform increase in $\tau$ with
increasing $\alpha$ for all $\kappa$ values (see Fig. 11(b)). The maximum
translocation times are observed when both $\kappa$ and $\alpha$ are large.

We look at the waiting time distribution $w(s)$ of the semiflexible polymer as
it moves from the cis to the trans side of the channel at low forces
($f_{0}=0.2$). In Fig. 12, we have plotted $w(s)$ for different values of the
bending rigidity as the apex angle is varied. For low bending rigidity, $w(s)$
shows an increase with $s$ initially and then a largely flat region at
intermediate $s$, for all values of $\alpha$ (Fig. 12(a)). These results are
similar to that for uniformly extended channels. Near the end, the attractive
interactions prevent exit of the polymer leading to larger waiting times. For
stiffer polymers, the waiting times for all monomers monotonically increase
with $\alpha$ which is consistent with the behavior of total translocation
time. As the polymer rigidity increases (Fig. 12(b)), the entropic gain in
moving towards the trans end of the conical channel becomes lesser. The
translocation dynamics of the end monomers is now dominated by the attractive
surface interactions. The end segment of the polymer keeps deflecting between
the walls and are stuck in either wall for long periods, thereby increasing
the waiting times.

##  IV Free Energy

(a)  (b)

(c)  (d)

Figure 13: (a) Entropic part of the free energy
($\mathscr{F}_{c}+\mathscr{F}_{e}$) for a flexible polymer as a function of
the translocation coordinate $s$ for different $\alpha$ values showing non-
monotonic behavior. (b) Energetic contribution to free energy for a flexible
polymer from surface interactions ($\mathscr{F}_{lj}$) as a function of
translocation coordinate $s$ for various $\alpha$ values. (c) Entropic part of
the free energy ($\mathscr{F}_{c}+\mathscr{F}_{e}$) for a semiflexible polymer
($\kappa=8$) as a function of the translocation coordinate $s$ for different
$\alpha$ values showing largely monotonic behavior. (Inset) The barrier height
decreases monotonically with $\alpha$. (d) Energetic contribution to free
energy for the semiflexible polymer from surface interactions
($\mathscr{F}_{lj}$) as a function of translocation coordinate $s$ for various
$\alpha$ values.

It has been established that the relaxation time of a polyelectrolyte chain in
equilibrium is much shorter than the translocation time for the typical range
of chain lengths used and typical voltages applied in most experiments Katkar
and Muthukumar (2018). This implies that a polymer chain initially at
equlibrium can be taken to be in quasi-equilibrium as it translocates across
the nanochannel especially for low external forces. In fact the quasi-
equlibrium approximation in such cases can be even extended to translocation
of every segment of the polymer Katkar and Muthukumar (2018). The quasi-
equlibrium approximation obviously fails for high forces and longer lengths of
the polymer. Using this approximation for low forces and chain lengths
considered in our simulations, we construct the free energy landscape for the
translocation process to understand the translocation behavior obtained in our
simulations.

There are several contributions to the free energy : (i) the entropic
contribution coming from the confinement of the polymer segment inside the
channel (ii) entropic contribution due to the polymer segments outside the
channel (ii) energetic contribution due to surface interactions between the
channel and the polymer and (iv) energetic contribution due to the external
force. The total free energy of the system in terms of the translocation
coordinate can be written as

$\mathscr{F}_{t}(s)=\mathscr{F}_{e}+\mathscr{F}_{c}+\mathscr{F}_{pp}$ (7)

where $\mathscr{F}_{c}$ is the confinement free energy, $\mathscr{F}_{e}$ is
the entropic component due to polymer segments outside the channel and
$\mathscr{F}_{pp}$ is the free energy contribution due to the channel polymer
interactions.

To construct $\mathscr{F}_{e}$ for a semiflexible polymer segment inside a
conical channel, we use the approach in Nikoofard  _et al._ (2013). The
contribution to the free-energy from a polymer of size $R$ confined in an area
$A$ of the cone in two-dimensions is given as $\mathscr{F}_{c}\sim
k_{B}T\left(\frac{R^{2}}{A}\right)^{\frac{1}{2\nu-1}}$ Sakaue and Raphaël
(2006). If $N_{p}$ denotes the number of polymer beads confined inside the
channel, then $R\sim\sigma N_{p}^{\nu}$, where $\nu=3/4$ is the Flory exponent
in two-dimensions. $A$ is dependent on the degree of confinement. If $L$
denotes the length of the region inside the channel along the $x-$direction
that holds the confined polymer and $a$ the length which is empty, then $A$ is
given as
$A=\frac{1}{4\tan\alpha}\left[(D_{0}+2(a+L)\tan\alpha)^{2}-(D_{0}+2a\tan\alpha)^{2}\right]$.

In the case of rigid polymers, the analysis can be extended by replacing the
number of monomers $N_{p}$ with the number of Kuhn segments
$\tilde{N_{p}}\approx N_{p}\sigma/l_{K}$ where $l_{K}=2l_{p}$ is the Kuhn
length. $l_{p}$ is the persistence length of the polymer. In two-dimensions,
the bending rigidity is related to the persistence length as
$\kappa/k_{B}T=l_{p}/2$. We would however like to express $\mathscr{F}_{e}$ in
terms of the translocation coordinate $s$. Therefore, we extract the number of
beads that are inside the channel from the simulations for a fixed
translocation coordinate. As is expected, there would be several
configurations of the polymer (with different values of $N_{p}$) inside the
channel for a given $s$. The confinement free energies of all such
configurations are averaged to give $\mathscr{F}_{c}(s)$.

The free energy contribution from the polymer segments outside the channel is
given asGauthier and Slater (2008b):

$\mathscr{F}_{e}=\begin{cases}\gamma^{\prime}\log N_{c}&\text{ $N_{c}>1$ and
$N_{t}=0$}\\\ \gamma^{\prime}\log N_{c}+\gamma^{\prime}\log N_{t}&\text{
$N_{c}>1$ and $N_{t}>1$}\\\ \gamma^{\prime}\log N_{t}&\text{ $N_{c}=0$ and
$N_{t}>1$}\end{cases}$

where $N_{c}$ is the number of beads on the cis-side of the channel and
$N_{T}$ is the number of beads on the trans-side of the channel. Note that
$N_{c}$, $N_{T}$ are functions of the translocation coordinate $s$ and are
obtained from the simulation similar to the confinement free energy for every
value of $s$. Further, $N_{c}+N_{p}+N_{T}=N$ at all times.

The free energy contribution from the surface interaction $\mathscr{F}_{lj}$
is obtained by summing over the LJ potential felt by each polymer bead when
they are inside the channel. For a given translocation coordinate, we identify
the number of beads inside the channel and calculate the total potential felt
by the beads. As we mentioned before, for a given $s$, there would be multiple
configurations of the polymer inside the channel (corresponding to the time
steps for which $s=1,2,...$). All such contributions are appropriately
accounted for and the averaged free energy contribution $\mathscr{F}_{lj}$ due
to surface interactions is evaluated numerically as the translocation
coordinate changes.

In Fig. 13(a), we plot the total free energy contributions which are entropic
in origin for a flexible polymer, arising due to the confined segments of the
polymer and the segments outside the channel. The free energy contribution due
to channel-polymer interactions are plotted in Fig. 13(b). The entropic
contribution exhibits a free energy barrier for an intermediate $s$ for all
values of the channel angle $\alpha$. This indicates that entropically it is
more difficult for the polymer to enter the channel from the cis side.
However, the energetic contribution coming from the attractive $LJ$ potential
dominates the initial stages and the polymer is sucked into the channel. At
intermediate $s$, the competition between the two leads to the comparatively
larger waiting times observed in Fig. 8. As $\alpha$ is increased, we observe
a clear non-monotonic behavior in the entropic part of the free energy. The
free energy barrier increases as we go from $\alpha=1.0^{\circ}$ to
$\alpha=2.7^{\circ}$ and decreases as $\alpha$ is increased further. This
indicates that translocation would be more difficult from $\alpha=1.0^{\circ}$
to $\alpha=2.7^{\circ}$ and favorable as the channel angle increases. This is
consistent with the behavior observed in the translocation time as a function
of the channel angle, where for $f_{0}=0.1$, $\tau$ shows a peak at
$\alpha\approx 2.7^{\circ}$ and decreases with increasing $\alpha$. The
energetic contribution $\mathscr{F}_{lj}$ does not show non-monotonic behavior
at these values of $\alpha$, with significant overlap in the free energy
profiles (see Fig. 13(b)). As $\alpha$ is increased, $\mathscr{F}_{lj}$
becomes shallower indicating easier translocation and therefore decreasing
$\tau$.

In Fig. 13(c), we plot the total free energy contributions which are entropic
in origin for a semiflexible polymer for various channel angles. Similar to
the flexible polymer, the entropic contribution indicates a barrier at
intermediate $s$. However, unlike the flexible polymer, this barrier reduces
monotonically with increasing $\alpha$. The corresponding free energy
contribution due to channel-polymer interactions also becomes shallower
monotonically as $\alpha$ is increased (see Fig. 13(d)). Therefore, for a
semiflexible polymer, translocation becomes easier with increasing $\alpha$.
If we compare our analytical predictions with that of the simulation results
in Fig. 9(a) for $\kappa=8$, we find that translocation time indeed decreases
with increasing $\alpha$ upto $\alpha=6^{\circ}$ consistent with our
predictions. However, for larger channel angles, the translocation time
increases. This behavior is very different from the predictions from the free
energy calculations.

In order to understand this behavior, it is crucial to look at the potential
landscape inside the channel. At large channel angles, the strength of the
potential due to the channel walls becomes negligible for large regions inside
the channel. Even during entry of the polymer from the cis side of the
channel, the attractive potential becomes much weaker in the region close to
the cis side for large angles. This hampers the pulling ability of the channel
resulting in larger waiting times of even the initial beads of the polymer. As
the polymer is pulled in, the segments of the stiff polymer inside the channel
keeps getting deflected from one wall to the other (see Supplementary Movie).
The potential near the centre of the channel is negligible as is the external
force and the polymer segments stay attached to either wall for considerable
times. This leads to increasing translocation times at higher channel angles.
For smaller channel angles, the average positions of the polymer segments are
much closer to the middle of the channel as both walls try to pull on them.
With external force much stronger near the center, the polymer is pulled out
faster. Note that for flexible polymers, the possibility of more conformations
of the polymer as compared to stiffer chains at larger angles ensures that the
situation is much more homogeneous and we do not observe much change at larger
channel angles.

##  V Discussions

The translocation of polymers through asymmetric shaped channels has sparked a
lot of interest due to its significant advantages in bio-sensing applications.
Buoyed by experimental studies using conical nanochannels, we have studied the
translocation dynamics of a semiflexible polymer through such a channel under
the influence of a spatially varying external drive and attractive surface
interactions. The waiting time distribution show rich features arising due to
the polymer stiffness, surface interactions and nature of the external drive.
We attempt to understand some of these features using a free energy argument
based on a quasi-equilibrium approximation which is applicable for smaller
polymer lengths and low forces. The theory captures some of the non-monotonic
features of the translocation dynamics. At higher stiffness and larger channel
angles, the variable potential landscape arising from the surface interactions
becomes dominant. This variability needs to be accounted for in order to
provide a more accurate description of the dynamics.

In our analysis we have ignored the role of electrostatic interactions
Buyukdagli  _et al._ (2019). Translocation rate depends on the charge density
of the polymer and the density of the bulk concentration of the solvent in a
confined pore. Low polymer charge density and low solvent concentration leads
to lower translocation rates, while translocation rate increases sharply for
high polymer charge density at high solvent concentration. It will be
interesting to ask how these rates could be modified when considering
asymmetric channels like the one considered in this study.

Recent experiments have revealed that driven polymer translocation through
synthetic nanopores is a two stage process with translocation initially
slowing with time before accelerating close to the end of the process Chen
_et al._ (2021). In our detailed simulations, we show however that the
translocation dynamics is strongly dependent on the stiffness of the polymer
and surface interactions. Attractive surface interactions can considerably
slow down the translocation process near the ends while increasing the angle
of the conical pore may facilitate the process. We have considered one
possible conical structure where the polymer enters the pore from the
constricted side of the channel. In future work we will explore the detailed
dynamics of the polymer as it enters from the wider end of the channel which
is expected to give rise to significantly different dynamics due to the
asymmetry of the channel.

##  VI ACKNOWLEDGEMENTS:

AS acknowledges University Grants Commission (UGC), India for the financial
support under the UGC-JRF/SRF scheme. We are thankful to the computational
facility provided by IISER, Mohali.

## References

  * Muthukumar (2016) M. Muthukumar,  _Polymer translocation_ (CRC press, 2016).
  * Palyulin  _et al._ (2014) V. V. Palyulin, T. Ala-Nissila, and R. Metzler, Soft matter 10, 9016 (2014).
  * Branton  _et al._ (2010) D. Branton, D. W. Deamer, A. Marziali, H. Bayley, S. A. Benner, T. Butler, M. Di Ventra, S. Garaj, A. Hibbs, X. Huang, _et al._ , in  _Nanoscience and technology: A collection of reviews from Nature Journals_ (World Scientific, 2010) pp. 261–268.
  * Wanunu (2012) M. Wanunu, Physics of life reviews 9, 125 (2012).
  * Howorka and Siwy (2009) S. Howorka and Z. Siwy, Chemical Society Reviews 38, 2360 (2009).
  * Keyser (2011) U. F. Keyser, Journal of The Royal Society Interface 8, 1369 (2011).
  * Milchev (2011) A. Milchev, Journal of Physics: Condensed Matter 23, 103101 (2011).
  * Bezrukov  _et al._ (1994) S. M. Bezrukov, I. Vodyanoy, and V. A. Parsegian, Nature 370, 279 (1994).
  * Kasianowicz  _et al._ (1996) J. J. Kasianowicz, E. Brandin, D. Branton, and D. W. Deamer, Proceedings of the National Academy of Sciences 93, 13770 (1996).
  * Sung and Park (1996) W. Sung and P. Park, Physical review letters 77, 783 (1996).
  * Park and Sung (1998) P. J. Park and W. Sung, The Journal of chemical physics 108, 3013 (1998).
  * Muthukumar (1999) M. Muthukumar, The Journal of Chemical Physics 111, 10371 (1999).
  * Muthukumar (2003) M. Muthukumar, The Journal of chemical physics 118, 5174 (2003).
  * Sakaue (2007) T. Sakaue, Physical Review E 76, 021803 (2007).
  * Sakaue (2010) T. Sakaue, Physical Review E 81, 041808 (2010).
  * Saito and Sakaue (2012) T. Saito and T. Sakaue, Physical Review E 85, 061803 (2012).
  * Sarabadani  _et al._ (2014) J. Sarabadani, T. Ikonen, and T. Ala-Nissila, The Journal of Chemical Physics 141, 214907 (2014).
  * Rowghanian and Grosberg (2011) P. Rowghanian and A. Y. Grosberg, The Journal of Physical Chemistry B 115, 14127 (2011).
  * Ikonen  _et al._ (2012a) T. Ikonen, A. Bhattacharya, T. Ala-Nissila, and W. Sung, Physical Review E 85, 051803 (2012a).
  * Ikonen  _et al._ (2013) T. Ikonen, A. Bhattacharya, T. Ala-Nissila, and W. Sung, EPL (Europhysics Letters) 103, 38001 (2013).
  * Ikonen  _et al._ (2012b) T. Ikonen, A. Bhattacharya, T. Ala-Nissila, and W. Sung, The Journal of chemical physics 137, 085101 (2012b).
  * Bhattacharya and Binder (2010) A. Bhattacharya and K. Binder, Physical Review E 81, 041804 (2010).
  * Bhattacharya (2013) A. Bhattacharya, Polymer Science Series C 55, 60 (2013).
  * de Haan and Slater (2010) H. W. de Haan and G. W. Slater, Physical Review E 81, 051802 (2010).
  * Gauthier and Slater (2008a) M. G. Gauthier and G. W. Slater, The Journal of chemical physics 128, 02B612 (2008a).
  * Gauthier and Slater (2008b) M. G. Gauthier and G. W. Slater, The Journal of Chemical Physics 128, 05B619 (2008b).
  * Huang  _et al._ (2014) A. Huang, A. Bhattacharya, and K. Binder, The Journal of Chemical Physics 140, 214902 (2014).
  * Luo  _et al._ (2007) K. Luo, T. Ala-Nissila, S.-C. Ying, and A. Bhattacharya, Physical Review Letters 99, 148102 (2007).
  * Luo  _et al._ (2008) K. Luo, T. Ala-Nissila, S.-C. Ying, and A. Bhattacharya, Physical Review E 78, 061918 (2008).
  * Cohen  _et al._ (2011) J. A. Cohen, A. Chaudhuri, and R. Golestanian, Physical review letters 107, 238102 (2011).
  * Cohen  _et al._ (2012a) J. A. Cohen, A. Chaudhuri, and R. Golestanian, The Journal of chemical physics 137, 204911 (2012a).
  * Sarabadani and Ala-Nissila (2018) J. Sarabadani and T. Ala-Nissila, Journal of Physics: Condensed Matter 30, 274002 (2018).
  * Suhonen and Linna (2018) P. M. Suhonen and R. P. Linna, Physical Review E 97, 062413 (2018).
  * Katkar and Muthukumar (2018) H. H. Katkar and M. Muthukumar, The Journal of chemical physics 148, 024903 (2018).
  * Kumar Sharma  _et al._ (2019) R. Kumar Sharma, I. Agrawal, L. Dai, P. S. Doyle, and S. Garaj, Nature communications 10, 1 (2019).
  * Suma and Micheletti (2017) A. Suma and C. Micheletti, Proceedings of the National Academy of Sciences 114, E2991 (2017).
  * Suma  _et al._ (2020) A. Suma, L. Coronel, G. Bussi, and C. Micheletti, Nature communications 11, 1 (2020).
  * Plesa  _et al._ (2016) C. Plesa, D. Verschueren, S. Pud, J. Van Der Torre, J. W. Ruitenberg, M. J. Witteveen, M. P. Jonsson, A. Y. Grosberg, Y. Rabin, and C. Dekker, Nature nanotechnology 11, 1093 (2016).
  * Schneider and Dekker (2012) G. F. Schneider and C. Dekker, Nature biotechnology 30, 326 (2012).
  * Sun and Luo (2014) L.-Z. Sun and M.-B. Luo, Journal of Physics: Condensed Matter 26, 415101 (2014).
  * Sun  _et al._ (2018) L.-Z. Sun, H. Li, X. Xu, and M.-B. Luo, Journal of Physics: Condensed Matter 30, 495101 (2018).
  * Wong and Muthukumar (2010) C. T. A. Wong and M. Muthukumar, The Journal of chemical physics 133, 07B607 (2010).
  * Meller  _et al._ (2001) A. Meller, L. Nivon, and D. Branton, Physical Review Letters 86, 3435 (2001).
  * Meller (2003) A. Meller, Journal of physics: condensed matter 15, R581 (2003).
  * Cohen  _et al._ (2012b) J. A. Cohen, A. Chaudhuri, and R. Golestanian, Physical Review X 2, 021002 (2012b).
  * Kumar  _et al._ (2018) R. Kumar, A. Chaudhuri, and R. Kapri, The Journal of chemical physics 148, 164901 (2018).
  * Derrington  _et al._ (2010) I. M. Derrington, T. Z. Butler, M. D. Collins, E. Manrao, M. Pavlenok, M. Niederweis, and J. H. Gundlach, Proceedings of the National Academy of Sciences 107, 16060 (2010).
  * Venkatesan and Bashir (2011) B. M. Venkatesan and R. Bashir, Nature nanotechnology 6, 615 (2011).
  * Dekker (2010) C. Dekker, Nanoscience And Technology: A Collection of Reviews from Nature Journals , 60 (2010).
  * Wanunu and Meller (2007) M. Wanunu and A. Meller, Nano letters 7, 1580 (2007).
  * Gershow and Golovchenko (2007) M. Gershow and J. A. Golovchenko, Nature nanotechnology 2, 775 (2007).
  * Luan  _et al._ (2012) B. Luan, G. Stolovitzky, and G. Martyna, Nanoscale 4, 1068 (2012).
  * Keyser  _et al._ (2006) U. F. Keyser, B. N. Koeleman, S. Van Dorp, D. Krapf, R. M. Smeets, S. G. Lemay, N. H. Dekker, and C. Dekker, Nature Physics 2, 473 (2006).
  * Mathé  _et al._ (2005) J. Mathé, A. Aksimentiev, D. R. Nelson, K. Schulten, and A. Meller, Proceedings of the National Academy of Sciences 102, 12377 (2005).
  * Harrell  _et al._ (2006) C. C. Harrell, Y. Choi, L. P. Horne, L. A. Baker, Z. S. Siwy, and C. R. Martin, Langmuir 22, 10837 (2006).
  * Lan  _et al._ (2011) W.-J. Lan, D. A. Holden, and H. S. White, Journal of the American Chemical Society 133, 13300 (2011).
  * Thacker  _et al._ (2012) V. V. Thacker, S. Ghosal, S. Hernández-Ainsa, N. A. Bell, and U. F. Keyser, Applied physics letters 101, 223704 (2012).
  * Zhou  _et al._ (2017) J. Zhou, Y. Wang, L. D. Menard, S. Panyukov, M. Rubinstein, and J. M. Ramsey, Nature communications 8, 1 (2017).
  * Chen  _et al._ (2021) K. Chen, I. Jou, N. Ermann, M. Muthukumar, U. F. Keyser, and N. A. Bell, Nature Physics 17, 1043 (2021).
  * Nikoofard  _et al._ (2013) N. Nikoofard, H. Khalilian, and H. Fazli, The Journal of Chemical Physics 139, 074901 (2013).
  * Nikoofard and Fazli (2015) N. Nikoofard and H. Fazli, Soft Matter 11, 4879 (2015).
  * Tu  _et al._ (2018) B. Tu, S. Bai, B. Lu, and Q. Fang, Scientific reports 8, 1 (2018).
  * Sun  _et al._ (2021) L.-Z. Sun, W.-P. Cao, C.-H. Wang, and X. Xu, The Journal of Chemical Physics 154, 054903 (2021).
  * Adhikari and Bhattacharya (2013) R. Adhikari and A. Bhattacharya, The Journal of Chemical Physics 138, 204909 (2013).
  * Sarabadani  _et al._ (2017) J. Sarabadani, T. Ikonen, H. Mökkönen, T. Ala-Nissila, S. Carson, and M. Wanunu, Scientific reports 7, 1 (2017).
  * Odijk (1983) T. Odijk, Macromolecules 16, 1340 (1983).
  * Sakaue and Raphaël (2006) T. Sakaue and E. Raphaël, Macromolecules 39, 2621 (2006).
  * Buyukdagli  _et al._ (2019) S. Buyukdagli, J. Sarabadani, and T. Ala-Nissila, Polymers 11, 118 (2019).


# Pedagogical Alignment of Large Language Models

Shashank Sonkar
Rice University
Houston, TX
<EMAIL_ADDRESS>
&Kangqi Ni 11footnotemark: 1
Rice University
Houston, TX
<EMAIL_ADDRESS>
Sapana Chaudhary
Texas A&M University
College Station, TX
<EMAIL_ADDRESS>
&Richard G. Baraniuk
Rice University
Houston, TX
<EMAIL_ADDRESS>
Equal contribution.

###### Abstract

In this paper, we introduce the novel concept of pedagogically aligned Large
Language Models (LLMs) that signifies a transformative shift in the
application of LLMs within educational contexts. Rather than providing direct
responses to user queries, pedagogically-aligned LLMs function as scaffolding
tools, breaking complex problems into manageable subproblems and guiding
students towards the final answer through constructive feedback and hints. The
objective is to equip learners with problem-solving strategies that deepen
their understanding and internalization of the subject matter. Previous
research in this field has primarily applied the supervised finetuning
approach without framing the objective as an alignment problem, hence not
employing reinforcement learning through human feedback (RLHF) methods. This
study reinterprets the narrative by viewing the task through the lens of
alignment and demonstrates how RLHF methods emerge naturally as a superior
alternative for aligning LLM behaviour. Building on this perspective, we
propose a novel approach for constructing a reward dataset specifically
designed for the pedagogical alignment of LLMs. We apply three state-of-the-
art RLHF algorithms and find that they outperform SFT significantly. Our
qualitative analyses across model differences and hyperparameter sensitivity
further validate the superiority of RLHF over SFT. Also, our study sheds light
on the potential of online feedback for enhancing the performance of
pedagogically-aligned LLMs, thus providing valuable insights for the
advancement of these models in educational settings.

##  1 Introduction

Large Language Models (LLMs) are an integral part of the modern digital
landscape, and ensuring their alignment to be honest, harmless, and helpful -
the ‘3H’ principle - is critical to their effective usage [1]. However, the
needs of a student learning a subject or preparing for an exam differ
substantially from those of a casual internet user. While the latter might
seek immediate answers to their queries, the former benefits more from a
structured, guided approach that encourages active learning [2] and critical
thinking [3]. Therefore, if we aim to leverage LLMs to assist students’
learning, we must redefine the ‘helpfulness’ criterion in a way that
specifically supports this demographic [4].

Figure 1: The image depicts the comparison between a traditional Large
Language Model (LLM) interaction (left) and a pedagogically-aligned LLM
interaction (right). The traditional LLM directly provides the user with the
answer, while the pedagogically-aligned LLM guides the student to the solution
by presenting a series of subproblems. This elucidates the concept of
pedagogical alignment, emphasizing the transformation from direct problem-
solving to a guided, scaffolded learning experience.

In this paper, we propose a novel concept of pedagogical alignment for LLMs
which aligns LLMs to function as scaffolding tools [5], breaking complex
problems into manageable subproblems and guiding students towards the final
answer through constructive feedback and hint [6]. This approach redefines the
‘helpfulness’ criterion of the ‘3H’ principle, tailoring it to the unique
needs and learning objectives of students. We argue that to truly assist
students, LLMs should not merely provide immediate answers, but rather guide
learners through a structured, step-by-step process that encourages active
learning and critical thinking. The pedagogical alignment of LLMs is not just
a theoretical construct; it is grounded in the belief that learning can be
significantly improved when students are given the opportunity to actively
participate in their learning process [7] guided by feedback [8].

Our work builds upon the framework of Conversational Learning with Analytical
Step-by-Step Strategies (CLASS) [9, 10], a tool that applies above principles
of learning sciences to construct an Intelligent Tutoring System (ITS) using
LLMs. However, their approach is implemented through supervised fine-tuning of
the LLM, which we identify as a limitation. We argue that the desired behavior
to foster learning can be more effectively modeled using reinforcement
learning. However, the challenge lies in the absence of a reward dataset
necessary for reinforcement learning in education context.

To address this challenge, we reinterpret the CLASS framework from a
pedagogical alignment perspective, proposing that it can be viewed as a
potential reinforcement learning problem where the rewards are associated with
pedagogical alignment. This perspective allows us to construct a preference
(reward) dataset, laying the groundwork to enhance the alignment of LLMs to
better serve educational purposes.

Next, using this reward dataset, we study the effective of three RLHF
algorithms: Direct Preference Optimization (DPO) [11], Identity Preference
Optimization (IPO) [12], and Kahneman Tversky Optimization (KTO) [13] for
pedagogical alignment of two state-of-the-art open-source LLMs - Mistral [14]
and Zephyr [15]. We measure the the efficacy of these algorithms on
pedagogical alignment by an accuracy metric. This metric quantifies the
alignment by considering whether the correct pedagogical action was taken by
the model - a crucial aspect of our pedagogical alignment approach.

Our findings reveal that DPO and KTO significantly enhance the performance of
both Mistral and Zephyr models, improving our pedagogical alignment accuracy
by 37 and 23 points respectively, compared to CLASS’ SFT approach. Also, we
examine the influence of the beta hyperparameter on the performance of these
alignment algorithms since beta is a critical parameter that significantly
affects the outcome of the alignment process.

Our contributions in this paper are threefold:

  1. 1.

We introduce the concept of “pedagogical alignment,” a new approach to align
LLMs in a manner that fosters optimal student learning.

  2. 2.

We reinterpret the CLASS framework from a pedagogical alignment perspective,
which allows us to create a reward dataset.

  3. 3.

Using the reward dataset, we demonstrate the effectiveness of reinforcement
learning-based alignment algorithms on state-of-the-art LLMs.

Through this work, we aim to contribute to the effective use of LLMs in
education, highlighting their potential as productive learning aids and
enhancing the learning experience for students.

Table 1: Inference prompt of CLASS framework.

⬇

Prompt: Instructions to act as a Tutorbot, an AI-powered chatbot designed to
help student with a question. Your goal is to break problem into subproblems
for the student, and walk the student to the final answer using the
subproblems and hints.

For each student response, think which category your response falls in, and
then frame the reply using this json structure:

{ "Evaluation": "a,b,c,d,e,f,g",

"Action Based on Evaluation": "1,2,3,4,5,6,..,12",

"Subproblem State": "x,y,z",

"Subproblem": "..",

"Tutorbot": "..", }

Evaluation of Student Response:

a) Correct Response

b) Incorrect Response

c) Partially Correct Response

d) Ambiguous or Unclear Response

e) ...

Tutorbot Actions:

If "a" is the evaluation, then:

1) Confirm the correct answer. Check for completeness of the answer. If
solution is incomplete, notify the student to complete the solution.

⬇

Prompt continues:

If "b" is the evaluation, then:

2) Notify the student about the mistake, Provide constructive feedback to
pinpoint the errors, Offer hints

3) Step in to provide a solution if the student is unable to answer even after
multiple attempts.

If "c" is the evaluation, then:

4) Acknowledge the accurate parts, Promptly notify the student about the
mistake, Provide constructive feedback to pinpoint the errors, Offer helpful
hints

5) Step in to provide a solution if the student is unable to answer after
multiple attempts

If "d" is the evaluation, then:

6) Actively seek clarification through relevant follow-up questions. Request
student to provide specific information

Subproblem State:

x) One of the subproblems is currently being solved

y) Subproblem finished, move to next subproblem

z) Subproblem finished, no next subproblem, problem finished

##  2 CLASS through the lens of pedagogical alignment

This section presents the concept of pedagogical alignment and uses it to
examine the CLASS framework. This examination highlights the potential of
application of reinforcement learning and paves the way for constructing a
reward dataset to enhance the pedagogical alignment of LLMs.

###  2.1 Pedagogical alignment

In the landscape of LLMs, the term “alignment" is a concept used to describe
the process of adapting LLMs to align with human preferences. In the context
of online conversational agents, alignment specifically refers to making an
LLM honest, harmless, and helpful to a typical internet user [1]. However, the
needs of a student are distinct from those of a regular internet user. As
illustrated in Figure 1, when an internet user poses a question to a LLM, the
LLM’s response should directly address the user’s query. This is the primary
use case for out-of-the-box LLMs like Mistral [14] and Llama [16], which are
designed for such interactions. However, to promote optimal student learning,
LLMs should guide the student to reach the correct answer through their own
reasoning by employing effective scaffolding [17], involving problem
decomposition, constructive feedback, and timely hints.

Here, we introduce our unique perspective of viewing the CLASS framework in
the context of pedagogical alignment. The CLASS tutorbot differentiates itself
from conventional out-of-the-box LLMs through its approach to student queries.
Instead of directly providing answers, the CLASS tutorbot takes on a more
interactive and guided approach. It begins by breaking the problem into
smaller, manageable subproblems. Each subproblem presents a learning
opportunity for the student, guiding them towards the final solution in a
step-by-step manner.

The CLASS framework achieves effective scaffolding through a defined discrete
action and state space. For every student response, the CLASS tutorbot
determines whether the student’s response is correct, incorrect, or ambiguous.
For instance, if the response is incorrect, an out-of-the-box LLM’s
instinctual response would be to provide the correct answer to the user.
However, the CLASS tutorbot’s initial reaction is to offer constructive
feedback to the student and allow the student to arrive at the right answer
independently. If the student is unable to do so, then the CLASS’s eventual
response is to step in and provide the correct answer (refer to action 2 and 3
in inference prompt 1).

Similarly, if the student’s response is partially correct, the CLASS
tutorbot’s first reaction is to inform the student about the correct aspects
and then highlight the errors, encouraging the student to correct their answer
independently (action 4 in prompt 1), as opposed to providing the complete
answer, which is a typical reaction of an out-of-the-box LLM.

Pedagogical alignment: We define the term ‘pedagogical alignment’ to describe
this shift in behavior that CLASS introduces, in contrast to an out-of-the-box
LLM, to promote learning [18]. While out-of-the-box LLMs are designed to be
helpful to general users, pedagogically aligned LLMs are tailored to foster
optimal learning experiences for students.

###  2.2 Reinforcement learning algorithms for alignment

While the CLASS framework used SFT for pedagogical alignment, recent
developments in the field of Reinforcement Learning from Human Feedback (RLHF)
present promising alternatives. RLHF refines LLM’s behavior based on human-
derived preferences or feedback, promoting a more profound congruence with
human ethical standards and objectives. The seminal work by [19] introduced
the concept of leveraging human preferences to steer the fine-tuning process
of language models, thereby laying the foundational stone for RLHF. Subsequent
advancements [20, 21] have led to notable enhancements in the RLHF pipeline,
augmenting alignment efficiency and overall model performance.

The RLHF methodology unfolds in a structured three-step process:

  1. 1.

Supervised Fine-Tuning (SFT): Initially, the LLM undergoes SFT using a task-
specific dataset. This step adjusts the model’s initial parameters to align
with the downstream task’s objectives, providing a tailored foundation for
further refinement.

  2. 2.

Preference Modeling: Following SFT, the model generates multiple responses
under varied generation parameters. Human evaluators then rank these responses
based on their preferences. This ranking process informs the development of a
preference model, which learns to assign scores to LLM outputs, effectively
translating human judgments into quantifiable metrics.

  3. 3.

Online Reinforcement Learning: Utilizing the preference model, LLMs are
further aligned through online reinforcement learning algorithms such as
Proximal Policy Optimization (PPO) [22]. This phase iteratively refines the
model’s parameters to maximize the alignment with human feedback, enhancing
the model’s performance and ethical alignment.

More recently, Direct Preference Optimization (DPO) [11] has emerged as a
streamlined and robust advancement over RLHF, offering superior performance by
forgoing the necessity for a reward model. DPO enhances alignment efficiency
by directly optimizing a preference-based loss function, thereby simplifying
the implementation and the operational efficiency. DPO is designed to directly
leverage a dataset of preferences, represented as ${(x,y_{w},y_{l})}$, where
each tuple consists of a prompt $x$, a preferred response $y_{w}$, and a dis-
preferred response $y_{l}$. Let, $\pi_{\theta}$ be the LLM being finetuned,
and $\pi_{\text{ref}}$ be a reference LLM (generally an SFTed model). Then,
DPO’s optimization problem can be written as:

$L_{\text{DPO}}(\pi_{\theta};\pi_{\text{ref}})=-\mathbb{E}_{(x,y_{w},y_{l})\sim
D}\left[\log\sigma\left(\beta\log\frac{\pi_{\theta}(y_{w}|x)}{\pi_{\text{ref}}(y_{w}|x)}-\beta\log\frac{\pi_{\theta}(y_{l}|x)}{\pi_{\text{ref}}(y_{l}|x)}\right)\right],$
(1)

where $\sigma$ is the sigmoid function, and $\beta$ is an algorithmic
parameter. However, DPO has been observed to over-fit to the preference
dataset, and is believed to be sensitive to hyperparameter adjustments. To
prevent the over-fitting behaviour of DPO’s optimization algorithm, recent
proposals include Identity Preference Optimisation (IPO) [12] and Kahneman-
Tversky Optimization (KTO) [13].

All the three methodologies, DPO, IPO, and KTO, require a dataset composed of
preferences, denoted as ${(x,y_{w},y_{l})}$. The assembly of such datasets in
a real-world context is notably labor-intensive and financially burdensome.
However, in the subsequent section, we introduce our innovative
reward/preference data generation approach, which effectively circumvents
these challenges, paving the way for efficient and cost-effective pedagogical
alignment of LLMs.

##  3 Preference Data Generation for Pedagogical Alignment

Figure 2: This illustration shows the process of generating pedagogically-
aligned preference data using the CLASS framework. The GPT-student asks a
question, to which both the GPT-tutor and SFT-tutor respond. The key focus
here is on the divergence in the ‘Action Based on Evaluation’ between the two
tutors. In this example, the GPT-tutor’s response is deemed more pedagogically
aligned because it encourages the student to engage in critical thinking and
attempt the problem again, instead of directly providing the correct answer.
This action mismatch between the tutor responses allows us to construct a
preference dataset that distinguishes between the pedagogically preferred
(chosen) and less effective (rejected) responses. The preference dataset
serves as the foundation for applying reinforcement learning-based alignment
algorithms, with the goal of enhancing the pedagogical alignment of large
language models.

In this section, we describe our DPO, IPO, and KTO compatible data generation
pipeline. The objective is to generate a preference dataset $\mathcal{D}_{p}$
which necessitates distinguishing between pedagogically chosen and rejected
tutor responses based on their alignment with desired scaffolding strategies.
$\mathcal{D}_{p}$ will mimic a dataset structure that is compatible with DPO,
IPO, and KTO to optimally represent these preferences. This dataset has three
fundamental components: context, chosen, and rejected. In our case, the
context represents the accumulated conversational history between the student
and the tutor leading up to a particular interaction point, providing the
necessary backdrop for the tutor’s next dialogue turn. The chosen contains the
pedagogically preferred or optimal tutor responses to the context, while
rejected includes those tutor responses considered less effective or
misaligned with the scaffolding approach.

CLASS Conversational Dataset $\mathcal{D}_{c}$ and Agents: Our base dataset is
the CLASS conversational dataset represented by
$\mathcal{D}_{c}=\\{(\mathbf{x}_{i},\mathbf{y}_{i})|i=1,\ldots,N\\}$,
comprising $N$ student-tutor conversational turns. Each $\mathbf{x}$
represents a sequence of student utterances, and each corresponding
$\mathbf{y}$ represents the tutor response.

To generate the preference dataset $\mathcal{D}_{p}$ from $\mathcal{D}_{c}$,
we employ three autonomous agents: GPT-student ($\mathcal{A}_{\text{G}}$),
GPT-tutor ($\mathcal{A}_{\text{T}}$), and SFT-tutor
($\mathcal{A}_{\text{S}}$). The conversation starts with GPT-student
$\mathcal{A}_{\text{G}}$ posing a question ($Q$) sourced from a question bank
to which both GPT-tutor $\mathcal{A}_{\text{T}}$ and SFT-tutor
$\mathcal{A}_{\text{S}}$ provide responses in the following JSON structure:

$\mathbf{y}_{i}=\left\\{\begin{array}[]{l}\text{``Evaluation of Student
Response": ``a, b, c, \ldots",}\\\ \text{``Action Based on Evaluation": ``1,
2, 3, \ldots, 12",}\\\ \text{``Subproblem State": ``x, y, z",}\\\
\text{``Subproblem": ``...",}\\\ \text{``Tutorbot": ``..."
}\end{array}\right.$ (2)

Defining Context: We define the context $C_{t}$ for any given response at time
$t$ as all preceding conversational turns up to that point:

$C_{t}=\\{(\mathbf{x}_{i}^{(t)},\mathbf{y}_{i}^{(t-1)})|i<t\\},$ (3)

where $\mathbf{x}_{i}^{(t)}$ denotes the sequence of student response till
time $t$, and $\mathbf{y}_{i}^{(t-1)}$ denotes all tutor responses up to, but
not including, the response at time $t$. This construction provides a detailed
backdrop for each next tutor interaction, feeding into the decision mechanism
for chosen and rejected responses.

Preference Dataset $\mathcal{D}_{p}$ Construction: Given this understanding of
context, $\mathcal{D}_{p}$ is structured to encapsulate interactions where the
pedagogically preferred (chosen) response diverges from the alternative
(rejected) based on specific pedagogical criteria. This divergence is
identified by a function $f$ that assesses the pivotal attributes of each
response. Hence, $\mathcal{D}_{p}$ is defined as:

$\displaystyle\mathcal{D}_{p}=\big{\\{}$
$\displaystyle(C_{t},\mathbf{y}_{t}^{\text{Chosen}},\mathbf{y}_{t}^{\text{Rejected}})\,|\,\mathbf{y}_{t}^{\text{Chosen}}=\mathbf{y}_{t}^{\mathcal{A}_{T}},\mathbf{y}_{t}^{\text{Rejected}}=\mathbf{y}_{t}^{\mathcal{A}_{S}},$
$\displaystyle\text{and }f(\mathbf{y}_{t}^{\mathcal{A}_{T}})\neq
f(\mathbf{y}_{t}^{\mathcal{A}_{S}}),\forall t\in\\{1,\ldots,N\\}\big{\\}}$ (4)

where $C_{t}=x,y_{w}=y_{choosen},y_{l}=y_{rejected}$ from DPO equation (1)
from section 2.2, and the function $f$ is the key to distinguishing between
chosen and rejected tutor responses. $f$ indicates a discrepancy between
$\mathcal{A}_{\text{T}}$’s and $\mathcal{A}_{\text{S}}$’s responses if any of
the fields “Evaluation", “Action Based on Evaluation", or “Subproblem State"
diverge between $\mathcal{A}_{\text{T}}$ and $\mathcal{A}_{\text{S}}$, $f$
asserts $\mathcal{A}_{\text{T}}$’s response as preferable, under the
assumption that $\mathcal{A}_{\text{T}}$ offers a more accurate pedagogical
model, alongside the rejected response from $\mathcal{A}_{\text{S}}$. This
condition is captured through the function $f$—mapping responses to their
pedagogical attributes:

$f(\mathbf{y}_{i}^{\mathcal{A}})=\left(\mathrm{Eval}(\mathbf{y}_{i}^{\mathcal{A}}),\,\mathrm{ActionEval}(\mathbf{y}_{i}^{\mathcal{A}}),\,\mathrm{SubprobState}(\mathbf{y}_{i}^{\mathcal{A}})\right)$
(5)

where,

  * •

$\text{Eval}(\mathbf{y}_{i}^{\mathcal{A}})$ extracts insights on the tutor’s
grading of student’s solution or question understanding.

  * •

$\text{ActionEval}(\mathbf{y}_{i}^{\mathcal{A}})$ outlines the tutor’s
recommended actions based on the evaluation.

  * •

$\text{SubprobState}(\mathbf{y}_{i}^{\mathcal{A}})$ provides an understanding
of the tutor’s perception of the student’s current understanding or progress
on sub-problems within a larger problem context.

The evaluation of these attributes by $f$ facilitates an objective basis to
deem one response as chosen (preferable) and another as rejected (less aligned
with desired pedagogical outcomes) through a systematic assessment of their
pedagogical value and alignment.

##  4 Experiments

In this section, we provide a comprehensive overview of our dataset
construction, experimental design, and the subsequent findings derived from
the application of state-of-the-art RLHF algorithms to train pedagogically-
aligned LLMs.

###  4.1 Dataset and Evaluation

Our experimental design requires an extensive conversational dataset to train
and test SFT and RLHF models effectively. The CLASS framework [9] used GPT-4
to generate 610 conversations to train their SFT model. However, we need to
generate a larger volume of conversations to train RLHF models and conduct
comprehensive testing. Following the strategy outlined by CLASS, we set about
creating these additional conversations. We used the CLASS scaffolding prompts
to generate the problems, subproblems, and associated hints, which were based
on biology topics from the OpenStax college book [23], consistent with CLASS’
context. These problems seed the context of the conversations generated
through GPT-4 in the next step using the CLASS conversational prompt. In these
simulated student-tutor conversations, the tutor uses the subproblems and
hints to guide the student toward the final answer. Through this methodology,
we successfully generated an additional 1128 conversations, resulting in 1738
conversations. Each conversation, on average, comprises approximately 8.2
rounds and 1311.0 words. This comprehensive dataset provides a solid platform
for training and comparing SFT and RLHF models.

To create the training and testing set for SFT and RLHF, we randomly sampled
our conversational dataset into 3 partitions: the SFT training dataset $D_{c}$
(600 conversations), the RLHF training dataset $D_{p}$ (600 conversations),
and the test dataset $D_{t}$ (450 conversations). As each round within the
conversations is processed into a separate data point, $D_{c}$ contains 4942
data points, $D_{p}$ contains 4921 data points, and $D_{t}$ contains 3701 data
points.

For the evaluation, we focus on the three classification fields within the
model JSON responses described in section 3: Evaluation of Student Response,
Action Based on Evaluation, and Subproblem State. We employ accuracy in these
three fields to measure the pedagogical alignment of LLMs since these metrics
offer a comprehensive view of the model’s performance, evaluating the
correctness of the responses, the appropriateness of the actions, and the
state of problem-solving. Hence, LLMs that achieve higher accuracy in these
fields are considered to provide more pedagogically aligned student
assistance.

Table 2: SFT vs DPO, IPO, KTO. Each cell shows the average accuracy across
three classification fields - Evaluation of student response, Action based on
the evaluation, and Subproblem State, respectively. We observed that for both
Mistral (Mistral-7B-Instruct-v0.2) and Zephyr (zephyr-7b-beta), all three RL
alignment algorithms - DPO, KTO, and IPO - showed notable accuracy gains in
accuracy across all fields as compared to SFT. For both models, DPO and KTO
outperformed IPO. IPO performed best for Vicuna (vicuna-7b-v1.5) compared to
the other two.

Model | SFT | DPO | IPO | KTO
---|---|---|---|---
Mistral | 35.0 (34,34,37) | 74.0 (72,72,78) | 57.7 (56,56,61) | 74.7 (73,73,80)
Vicuna |  74.3 (73,70,80) |  75.0 (73,72,80) |  76.0 (74,72,82) |  75.3 (73,72,81)
Zephyr | 10.7 (12,11,9) | 74.0 (72,71,79) | 21.7 (22,22,21) | 74.7 (73,72,79)

Table 3: LLM Performance on key benchmarks. This table presents the
performance of Mistral (Mistral-7B-Instruct-v0.2), Zephyr (zephyr-7b-beta),
and Vicuna(vicuna-7b-v1.5) on different NLP benchmarks [24, 25]. We selected
Mistral and Zephyr for our experiments based on their superior performance
across various benchmarks, which provided a thorough and fair assessment
across a spectrum of reasoning and general knowledge tasks. We chose Vicuna
because it served as the base model for the CLASS framework.

Model | Average | ARC | HellaSwag | MMLU | TruthfulQA | Winogrande | GSM8K
---|---|---|---|---|---|---|---
Mistral | 65.71 | 63.14 | 84.88 | 60.78 | 68.26 | 77.19 | 40.03
Zephyr | 61.95 | 62.03 | 84.36 | 61.07 | 57.45 | 77.74 | 29.04
Vicuna | 52.06 | 53.24 | 77.39 | 50.82 | 50.33 | 72.06 | 8.11

###  4.2 Models and and Training Details

We experimented on 3 different LLMs: Mistral (Mistral-7B-Instruct-v0.2),
Zephyr (zephyr-7b-beta), and Vicuna (vicuna-7b-v1.5). We selected Mistral and
Zephyr based on the their top performance on diverse NLP benchmarks (Table 3).
Among all, several benchmark in particular evaluate scientific knowledge and
education-based question answering capabilities: ARC (multiple choice
questions from grade school science exams) [26], MMLU (multiple choice
questions from 57 high school and college subjects across STEM, the
humanities, the social sciences, and more) [27], and GSM8K (linguistically
diverse grade school math word problems) [28]. These benchmarks offer a
comprehensive, unbiased evaluation on the models’ performance, especially in
education setting. Furthermore, despite Vicuna’s inferior performance compared
to Mistral and Zephyr, we included it in our selection because it was the base
model for the CLASS framework’s initial experiments and represented the state
of the art at its launch.

Our experiments included training and evaluating SFT and RLHF models. It took
around 22 NVIDIA A6000 GPU hours for one cycle (SFT + Evaluation + RLHF +
Evaluation). A total of 17 cycles, which are 374 GPU hours, are needed to
complete all the experiments. We conducted the SFT experiments using FastChat
[29]. For the choice of SFT hyperparameters, we refered to parameters used to
instruct fine-tune Vicuna [29]: We employed a learning rate of 2e-5 with the
AdamW optimizer [30], a batch size of 16, a cosine scheduler, a weight decay
of 0.05, a warmup ratio of 0.1, and a total number of 3 epochs. We conducted
the RLHF experiments using TRL [31]. For the choice of RLHF hyperparamters, we
referred to parameters used to align Zephyr [32]: We employed a learning rate
of 1e-7 with the AdamW optimizer, a batch size of 16, a cosine scheduler, a
weight decay of 0.05, a warmup ratio of 0.1, and a total number of 3 epochs.

###  4.3 Experiment 1: SFT vs DPO, IPO, and KTO

In this experiment, we compared three alignment algorithms: Direct Preference
Optimization (DPO), Identity Preference Optimization (IPO), and Kahneman
Tversky Optimization (KTO) on the SFT models. We measured the accuracy across
three classification fields in the responses: Evaluation of Student Response,
Action based on Evaluation, and Subproblem State.

We used the same RLHF training setting as described in Section 4.2 with a beta
of 0.1 for DPO, IPO, and KTO across all models. After applying the alignment
algorithms, we observed significant accuracy improvements on Mistral and
Zephyr SFT models. Specifically, Mistral demonstrated improvements of 39%,
22.7%, and 39.7% from DPO, IPO, and KTO respectively, while Zephyr showed
improvements of 63.3%, 11%, and 64% respectively (Table 2). Notably, for both
Mistral and Zephyr, DPO and KTO significantly outperformed IPO. However, the
Vicuna model, already noted for its high SFT performance, showed only marginal
improvements from the alignment; all three algorithms exhibited roughly
similar performance. Our results are in alignment with the findings of similar
studies [32], which have reported superior performance of DPO and KTO over IPO
on Zephyr when the model was fine-tuned on preference datasets, such as
Intel’s ORCA dataset pairs [33] and the Ultra-Feedback dataset [15].

This experiment underscored the pivotal role of alignment algorithms in
improving the accuracy of SFT models. Mistral and Zephyr, the top performing
models, exhibited pronounced improvements over Vicuna. Therefore, we
demonstrated the efficacy of pedagogical alignment on state of the art models.

###  4.4 Experiment 2: Effect of Beta on DPO and IPO

Figure 3: DPO and IPO performance as a function of beta. This figure depicts
how the beta hyper-parameter affects the average classification accuracy of
DPO and IPO aligned models. We found that both DPO and IPO performances were
significantly influenced by beta adjustments. Specifically, a beta of 0.3
optimized results for DPO, with Mistral and Zephyr achieving around 75%
accuracy, while a higher beta of 0.9 notably reduced their accuracy. In
contrast, IPO performed best with a beta of 0.1, and increasing the beta
decreased the accuracy.

In this experiment, we analyzed the impact of the beta parameter on the
performance of DPO and IPO alignment algorithms. Beta is a key hyperparameter
in all three alignment algorithms (DPO, IPO, and KTO). It controls the
strength of the Kullback-Leibler (KL) divergence penalty between the trained
model and the reference model within the loss function. Nevertheless, the
choice of beta depends on the specific model and dataset and thus necessitates
an empirical study.

Similar to experiments in Section 4.3, we measured the average accuracy across
the three classification fields. We found that the performance of DPO and IPO
is sensitive to the choice of beta. As depicted in figure 3, for DPO, we
observed that a beta value of 0.3 yielded the best results: Mistral and Zephyr
achieved an average accuracy of 74.7% and 75%. Conversely, a larger beta
value, such as 0.9, diminished the performance of the models. Zephyr suffered
a large loss of 25.7 % in average accuracy from using 0.9 instead of 0.3. As
opposed to IPO, we observed that a beta of 0.1 produced the best performance
(figure 3). A larger choice of beta in IPO almost always came with a loss of
performance. As beta increases, Mistral’s accuracy dropped from 57.7% to 36 %,
while Zephyr’s accuracy dropped from 21.7 % to 11.3 %. Our results corroborate
the evidence from related studies [32], demonstrating that various alignment
algorithms exhibit distinct sensitivities to beta values, necessitating
tailored beta selections.

##  5 Conclusion

In this paper, we have introduced and explored the concept of pedagogical
alignment of LLMs to align these models with educational objectives that are
aimed at fostering optimal student learning. By reinterpreting the CLASS
framework as a reinforcement learning problem, we constructed a reward
dataset, thereby providing a foundation for enhancing the pedagogical
alignment of LLMs. Our experimental results, derived from applying
reinforcement learning-based alignment algorithms on state-of-the-art open-
source LLMs, corroborate our hypothesis. The RLHF methods we employed have
proven to outperform the SFT approach significantly, favoring the alignment of
LLMs with pedagogical objectives. Furthermore, our study highlights the
potential of online feedback for enhancing the performance of pedagogically-
aligned LLMs. In conclusion, while our research has made significant strides
in aligning LLMs with educational needs, further exploration and refinement of
reinforcement learning methods for this purpose remain necessary. We
anticipate that our work will serve as a stepping stone towards the optimal
integration of LLMs in educational settings, paving the way for future
research in this exciting domain.

## Acknowledgments

This work was supported by NSF grants 1842378, ONR grant N0014-20-1-2534,
AFOSR grant FA9550-22-1-0060, and a Vannevar Bush Faculty Fellowship, ONR
grant N00014-18-1-2047.

## References

  * [1] Amanda Askell, Yuntao Bai, Anna Chen, Dawn Drain, Deep Ganguli, Tom Henighan, Andy Jones, Nicholas Joseph, Ben Mann, Nova DasSarma, et al.  A general language assistant as a laboratory for alignment.  arXiv preprint arXiv:2112.00861, 2021.
  * [2] Michelene TH Chi.  Active-constructive-interactive: A conceptual framework for differentiating learning activities.  Topics in cognitive science, 1(1):73–105, 2009.
  * [3] Richard Paul and Linda Elder.  Critical Thinking: How to Prepare Students for a Rapidly Changing World.  Foundation for Critical Thinking, 1995.
  * [4] Richard E. Mayer.  Multimedia Learning.  Cambridge University Press, 2001.
  * [5] Richard E. Mayer.  Should there be a three-strikes rule against pure discovery learning?  American psychologist, 59(1):14, 2004.
  * [6] John D. Bransford, Ann L. Brown, and Rodney R. Cocking, editors.  How People Learn: Brain, Mind, Experience, and School: Expanded Edition.  National Academy Press, 2000.
  * [7] Michael Prince.  Does active learning work? a review of the research.  Journal of Engineering Education, 93(3):223–231, 2004.
  * [8] John Hattie and Helen Timperley.  The power of feedback.  Review of Educational Research, 77(1):81–112, 2007.
  * [9] Shashank Sonkar, Naiming Liu, Debshila Mallick, and Richard Baraniuk.  CLASS: A Design Framework for Building Intelligent Tutoring Systems Based on Learning Science principles.  In Findings of the Association for Computational Linguistics: EMNLP 2023, pages 1941–1961, 2023.
  * [10] Shashank Sonkar, MyCo Le, Xinghe Chen, Naiming Liu, Debshila Basu Mallick, and Richard G Baraniuk.  Code Soliloquies for Accurate Calculations in Large Language Models.  arXiv preprint arXiv:2309.12161, 2023.
  * [11] Rafael Rafailov, Archit Sharma, Eric Mitchell, Stefano Ermon, Christopher D Manning, and Chelsea Finn.  Direct preference optimization: Your language model is secretly a reward model.  arXiv preprint arXiv:2305.18290, 2023.
  * [12] Mohammad Gheshlaghi Azar, Mark Rowland, Bilal Piot, Daniel Guo, Daniele Calandriello, Michal Valko, and Rémi Munos.  A general theoretical paradigm to understand learning from human preferences.  arXiv preprint arXiv:2310.12036, 2023.
  * [13] Kawin Ethayarajh, Winnie Xu, Niklas Muennighoff, Dan Jurafsky, and Douwe Kiela.  Kto: Model alignment as prospect theoretic optimization, 2024.
  * [14] Albert Q. Jiang, Alexandre Sablayrolles, Arthur Mensch, Chris Bamford, Devendra Singh Chaplot, Diego de las Casas, Florian Bressand, Gianna Lengyel, Guillaume Lample, Lucile Saulnier, Lélio Renard Lavaud, Marie-Anne Lachaux, Pierre Stock, Teven Le Scao, Thibaut Lavril, Thomas Wang, Timothée Lacroix, and William El Sayed.  Mistral 7b, 2023.
  * [15] Lewis Tunstall, Edward Beeching, Nathan Lambert, Nazneen Rajani, Kashif Rasul, Younes Belkada, Shengyi Huang, Leandro von Werra, Clémentine Fourrier, Nathan Habib, Nathan Sarrazin, Omar Sanseviero, Alexander M. Rush, and Thomas Wolf.  Zephyr: Direct distillation of lm alignment, 2023.
  * [16] Hugo Touvron, Thibaut Lavril, Gautier Izacard, Xavier Martinet, Marie-Anne Lachaux, Timothée Lacroix, Baptiste Rozière, Naman Goyal, Eric Hambro, Faisal Azhar, et al.  Llama: Open and efficient foundation language models.  arXiv preprint arXiv:2302.13971, 2023.
  * [17] Nisan Stiennon, Long Ouyang, Jeffrey Wu, Daniel Ziegler, Ryan Lowe, Chelsea Voss, Alec Radford, Dario Amodei, and Paul F Christiano.  Learning to summarize with human feedback.  Advances in Neural Information Processing Systems, 33:3008–3021, 2020.
  * [18] C Addison Stone.  The metaphor of scaffolding: Its utility for the field of learning disabilities.  Journal of learning disabilities, 31(4):344–364, 1998.
  * [19] Daniel M Ziegler, Nisan Stiennon, Jeffrey Wu, Tom B Brown, Alec Radford, Dario Amodei, Paul Christiano, and Geoffrey Irving.  Fine-tuning language models from human preferences.  arXiv preprint arXiv:1909.08593, 2019.
  * [20] Reiichiro Nakano, Jacob Hilton, Suchir Balaji, Jeff Wu, Long Ouyang, Christina Kim, Christopher Hesse, Shantanu Jain, Vineet Kosaraju, William Saunders, et al.  Webgpt: Browser-assisted question-answering with human feedback.  arXiv preprint arXiv:2112.09332, 2021.
  * [21] Long Ouyang, Jeffrey Wu, Xu Jiang, Diogo Almeida, Carroll Wainwright, Pamela Mishkin, Chong Zhang, Sandhini Agarwal, Katarina Slama, Alex Ray, et al.  Training language models to follow instructions with human feedback.  Advances in Neural Information Processing Systems, 35:27730–27744, 2022.
  * [22] John Schulman, Filip Wolski, Prafulla Dhariwal, Alec Radford, and Oleg Klimov.  Proximal policy optimization algorithms.  arXiv preprint arXiv:1707.06347, 2017.
  * [23] Mary Ann Clark, Matthew Douglas, and Jung Choi.  Biology 2e.  OpenStax, 2018.
  * [24] Edward Beeching, Clémentine Fourrier, Nathan Habib, Sheon Han, Nathan Lambert, Nazneen Rajani, Omar Sanseviero, Lewis Tunstall, and Thomas Wolf.  Open llm leaderboard, 2023.
  * [25] Leo Gao, Jonathan Tow, Baber Abbasi, Stella Biderman, Sid Black, Anthony DiPofi, Charles Foster, Laurence Golding, Jeffrey Hsu, Alain Le Noac’h, Haonan Li, Kyle McDonell, Niklas Muennighoff, Chris Ociepa, Jason Phang, Laria Reynolds, Hailey Schoelkopf, Aviya Skowron, Lintang Sutawika, Eric Tang, Anish Thite, Ben Wang, Kevin Wang, and Andy Zou.  A framework for few-shot language model evaluation, 12 2023.
  * [26] Peter Clark, Isaac Cowhey, Oren Etzioni, Tushar Khot, Ashish Sabharwal, Carissa Schoenick, and Oyvind Tafjord.  Think you have solved question answering? try arc, the ai2 reasoning challenge, 2018.
  * [27] Dan Hendrycks, Collin Burns, Steven Basart, Andy Zou, Mantas Mazeika, Dawn Song, and Jacob Steinhardt.  Measuring massive multitask language understanding, 2021.
  * [28] Karl Cobbe, Vineet Kosaraju, Mohammad Bavarian, Mark Chen, Heewoo Jun, Lukasz Kaiser, Matthias Plappert, Jerry Tworek, Jacob Hilton, Reiichiro Nakano, Christopher Hesse, and John Schulman.  Training verifiers to solve math word problems, 2021.
  * [29] Lianmin Zheng, Wei-Lin Chiang, Ying Sheng, Siyuan Zhuang, Zhanghao Wu, Yonghao Zhuang, Zi Lin, Zhuohan Li, Dacheng Li, Eric Xing, et al.  Judging llm-as-a-judge with mt-bench and chatbot arena.  arXiv preprint arXiv:2306.05685, 2023.
  * [30] Ilya Loshchilov and Frank Hutter.  Decoupled weight decay regularization, 2019.
  * [31] Leandro von Werra, Younes Belkada, Lewis Tunstall, Edward Beeching, Tristan Thrush, Nathan Lambert, and Shengyi Huang.  TRL: Transformer Reinforcement Learning, 2020.
  * [32] Kashif Rasul, Edward Beeching, Lewis Tunstall, Leandro von Werra, and Omar Sanseviero.  Preference tuning llms with direct preference optimization methods, 2024.
  * [33] Wing Lian, Bleys Goodson, Eugene Pentland, Austin Cook, Chanvichet Vong, and "Teknium".  Openorca: An open dataset of gpt augmented flan reasoning traces, 2023.