$\mathbf{b}_{L}=E[\boldsymbol{1}\\{D>h_{N}\\}D^{-1}(\mathbf{X}^{\prime}\mathbf{X})^{\ast}\mathbf{X}^{\prime}(\mathbf{Y}-\mathbf{W}\boldsymbol{\delta})]+\boldsymbol{\gamma}\mathbf{h},$

where $\mathbf{h}=E[\boldsymbol{1}\\{D\leq
h_{N}\\}(\begin{array}[]{cccc}1&\cdots&D^{L-1}\end{array})^{\prime}]$,
$\mathbf{m}(u)=E[(\mathbf{X}^{\prime}\mathbf{X})^{\ast}\mathbf{X}^{\prime}(\mathbf{Y}-\mathbf{W}\boldsymbol{\delta})\mid
D=u]$, and
$\boldsymbol{\gamma}=(\begin{array}[]{cccc}\mathbf{m}^{(1)}(0)/1!&\cdots&\mathbf{m}^{(L)}(0)/L!\end{array})$.
The proof is the same as that of Theorem 1.

######  Theorem 3.

Under $T>p$ and Assumptions 1, 5, and 6,
$\mathbf{b}_{L}=\boldsymbol{\beta}+O\left(E[D^{L}\boldsymbol{1}\\{D\leq
h_{N}\\}]\right)$.

The second result is the asymptotic normality of $\boldsymbol{\hat{\theta}}$
as an estimator for $\boldsymbol{\theta}$.

######  Theorem 4.

Suppose there is some $c>0$ such that

$Nh_{N}^{2+c}\rightarrow\infty\mbox{ and }Nh_{N}^{2(L+1)}\rightarrow 0$ (24)

as $h_{N}\rightarrow 0$ and $N\rightarrow\infty$. Under $T>p$ and Assumptions
1, 5, 6, and 7,

$\sqrt{N}E_{N}\left[\boldsymbol{\hat{\zeta}}\boldsymbol{\hat{\zeta}}^{\prime}\right]^{-1/2}(\boldsymbol{\hat{\theta}}-\boldsymbol{\theta})\rightarrow_{d}\mathcal{N}\left(0,\mathbf{I}_{p}\right)$

as $h_{N}\rightarrow 0$ and $N\rightarrow\infty$.

A proof of this theorem is presented in the subsequent subsection.

###  C.2 Proof of Theorem 4

We use the following notations in the proof of Theorem 4. We use
$\pi_{h}=P(D\leq h_{N})$ and $\pi_{h\setminus 0}=P(0<D\leq h_{N})$. In the
following proof, we characterize the influence function as
$\boldsymbol{\zeta}=\boldsymbol{\xi}_{1}+\frac{1}{\pi_{h\setminus
0}}\boldsymbol{\xi}_{2}\boldsymbol{\Pi}^{-1}\boldsymbol{\Psi}_{1:L}^{-1}\mathbf{h}+\mathbf{Q}\mathbf{{V}}^{-1}\boldsymbol{\xi}_{3}$,
where
$\boldsymbol{\Psi}_{1:L}=\mathrm{diag}\left(\begin{array}[]{cccc}h_{N}&\cdots&h_{N}^{L}\end{array}\right)$,

$\displaystyle\boldsymbol{\Pi}=$ $\displaystyle\frac{1}{\pi_{h\setminus
0}}\boldsymbol{\Psi}_{1:L}^{-1}E\left[\mathbf{D}_{1:L}\mathbf{D}_{1:L}^{\prime}\right]\boldsymbol{\Psi}_{1:L}^{-1},$
$\displaystyle\mathbf{Q}=$
$\displaystyle\mathbf{R}-E\left[\left(\boldsymbol{1}\\{D>h_{N}\\}D^{-1}+\mathbf{D}_{1:L}^{\prime}E\left[\mathbf{D}_{1:L}\mathbf{D}_{1:L}^{\prime}\right]^{-1}\mathbf{h}\right)(\mathbf{X}^{\prime}\mathbf{X})^{\ast}\mathbf{X}^{\prime}\mathbf{W}\right],$
$\displaystyle\boldsymbol{\xi}_{1}=$
$\displaystyle\boldsymbol{1}\\{D>h_{N}\\}D^{-1}(\mathbf{X}^{\prime}\mathbf{X})^{\ast}\mathbf{X}^{\prime}(\mathbf{Y}-\mathbf{W}\boldsymbol{\delta})-E[\boldsymbol{1}\\{D>h_{N}\\}D^{-1}(\mathbf{X}^{\prime}\mathbf{X})^{\ast}\mathbf{X}^{\prime}(\mathbf{Y}-\mathbf{W}\boldsymbol{\delta})],$
$\displaystyle\boldsymbol{\xi}_{2}=$
$\displaystyle((\mathbf{X}^{\prime}\mathbf{X})^{\ast}\mathbf{X}^{\prime}(\mathbf{Y}-\mathbf{W}\boldsymbol{\delta})-\boldsymbol{\gamma}\mathbf{D}_{1:L})\mathbf{D}_{1:L}^{\prime}\boldsymbol{\Psi}_{1:L}^{-1},\mbox{
and }$ $\displaystyle\boldsymbol{\xi}_{3}=$
$\displaystyle\boldsymbol{1}\\{D>h_{N}\\}\mathbf{W}\mathbf{M}_{\mathbf{X}}(\mathbf{Y}-\mathbf{W}\boldsymbol{\delta}).$

Their sample counterparts are

$\displaystyle\boldsymbol{\hat{\Pi}}=$ $\displaystyle\frac{1}{\pi_{h\setminus
0}}\boldsymbol{\Psi}_{1:L}^{-1}E_{N}\left[\mathbf{D}_{1:L}\mathbf{D}_{1:L}^{\prime}\right]\boldsymbol{\Psi}_{1:L}^{-1},$
$\displaystyle\boldsymbol{\hat{\xi}}_{1}=$
$\displaystyle\boldsymbol{1}\\{D>h_{N}\\}D^{-1}(\mathbf{X}^{\prime}\mathbf{X})^{\ast}\mathbf{X}^{\prime}(\mathbf{Y}-\mathbf{W}\boldsymbol{\hat{\delta}})-E_{N}[\boldsymbol{1}\\{D>h_{N}\\}D^{-1}(\mathbf{X}^{\prime}\mathbf{X})^{\ast}\mathbf{X}^{\prime}(\mathbf{Y}-\mathbf{W}\boldsymbol{\hat{\delta}})],$
$\displaystyle\boldsymbol{\hat{\xi}}_{2}=$
$\displaystyle((\mathbf{X}^{\prime}\mathbf{X})^{\ast}\mathbf{X}^{\prime}(\mathbf{Y}-\mathbf{W}\boldsymbol{\hat{\delta}})-\boldsymbol{\hat{\gamma}}\mathbf{D}_{1:L})\mathbf{D}_{1:L}^{\prime}\boldsymbol{\Psi}_{1:L}^{-1},\mbox{
and }$ $\displaystyle\boldsymbol{\hat{\xi}}_{3}=$
$\displaystyle\boldsymbol{1}\\{D>h_{N}\\}\mathbf{W}\mathbf{M}_{\mathbf{X}}(\mathbf{Y}-\mathbf{W}\boldsymbol{\hat{\delta}}).$

Except for $\boldsymbol{\Psi}_{1:L}$, they are all different from the
corresponding ones in the proof of Theorem 2.

Unlike Section 4, we do not assume that $\phi$ is bounded away from zero here.
Thus, the asymptotic properties for $(\pi_{h\setminus 0},\pi_{h})$ depends on
the tail index $\alpha$. By Assumption 5 (iv),

$\alpha>1\implies\pi_{h\setminus 0}=\pi_{h}.$

By Assumption 5 (iii), the term $\pi_{h\setminus 0}/h_{N}^{\alpha}$ is bounded
away from zero and infinity.

Using the auxiliary results in Section D, we are now ready to prove Theorem 4.

###### Proof of Theorem 4.

First, we are going to show

$E_{N}\left[\boldsymbol{\hat{\zeta}}\boldsymbol{\hat{\zeta}}^{\prime}\right]^{-1}Var(\boldsymbol{\zeta})=\mathbf{I}_{p}+o_{p}(1).$
(25)

By (24) and Lemmas 23 and 25,

$\displaystyle
E_{N}[\|\boldsymbol{\hat{\zeta}}-\boldsymbol{\zeta}\|^{2}]^{1/2}=$
$\displaystyle o_{p}\left(1\\{\alpha\leq 2\\}{\frac{\pi_{h}}{\pi_{h\setminus
0}}}+1\\{\alpha>2\\}\right),$ $\displaystyle
E_{N}[\|\boldsymbol{\zeta}\|^{2}]^{1/2}=$ $\displaystyle O\left(1\\{\alpha\leq
2\\}{\frac{\pi_{h}}{h_{N}^{2}}}+1\\{\alpha>2\\}\right),$
$\displaystyle(E_{N}-E)\left[\boldsymbol{\zeta}\boldsymbol{\zeta}^{\prime}\right]=$
$\displaystyle o_{p}\left(1\\{\alpha\leq
2\\}{\frac{\pi_{h}^{2}}{h_{N}^{2}\pi_{h\setminus
0}}}+1\\{\alpha>2\\}\right),\qquad\text{and}$ $\displaystyle
E[\boldsymbol{\zeta}]=$ $\displaystyle o(1).$

Therefore,

$\displaystyle\|E_{N}\left[\boldsymbol{\hat{\zeta}}\boldsymbol{\hat{\zeta}}^{\prime}\right]-Var(\boldsymbol{\zeta})\|\leq$
$\displaystyle
E_{N}[\|\boldsymbol{\hat{\zeta}}-\boldsymbol{\zeta}\|^{2}]+2E_{N}[\|\boldsymbol{\hat{\zeta}}-\boldsymbol{\zeta}\|^{2}]^{1/2}E_{N}[\|\boldsymbol{\zeta}\|^{2}]^{1/2}+\|(E_{N}-E)\left[\boldsymbol{\zeta}\boldsymbol{\zeta}^{\prime}\right]\|+\|E[\boldsymbol{\zeta}]\|^{2}$
$\displaystyle=$ $\displaystyle o_{p}\left(1\\{\alpha\leq
2\\}{\frac{\pi_{h}^{2}}{h_{N}^{2}\pi_{h\setminus 0}}}+1\\{\alpha>2\\}\right).$

By Lemma 24,

$E_{N}\left[\boldsymbol{\hat{\zeta}}\boldsymbol{\hat{\zeta}}^{\prime}\right]^{-1}Var(\boldsymbol{\zeta})-\mathbf{I}_{p}=-\left(Var(\boldsymbol{\zeta})+\left(E_{N}\left[\boldsymbol{\hat{\zeta}}\boldsymbol{\hat{\zeta}}^{\prime}\right]-Var(\boldsymbol{\zeta})\right)\right)^{-1}\left(E_{N}\left[\boldsymbol{\hat{\zeta}}\boldsymbol{\hat{\zeta}}^{\prime}\right]-Var(\boldsymbol{\zeta})\right)=o_{p}(1).$

Second, we are going to show

$\sqrt{N}(\boldsymbol{\hat{\theta}}-\boldsymbol{\theta}-(E_{N}-E)[\boldsymbol{\zeta}])=o_{p}\left(1\\{\alpha\leq
2\\}\sqrt{\frac{\pi_{h}^{2}}{h_{N}^{2}\pi_{h\setminus
0}}}+1\\{\alpha>2\\}\right).$ (26)

Note that

$\displaystyle\|\boldsymbol{\hat{\theta}}-\boldsymbol{\theta}-(E_{N}-E)[\boldsymbol{\zeta}]\|\leq$
$\displaystyle\|\mathbf{b}_{L}-\boldsymbol{\beta}+E[\boldsymbol{\zeta}]\|$
$\displaystyle+\frac{1}{\pi_{h\setminus
0}}\left(\left\|\boldsymbol{\hat{\Pi}}^{-1}-\boldsymbol{\Pi}^{-1}\right\|\|\boldsymbol{\Psi}_{1:L}^{-1}\mathbf{\hat{h}}\|+\|\boldsymbol{\Pi}^{-1}\|\|\boldsymbol{\Psi}_{1:L}^{-1}\left(\mathbf{\hat{h}}-\mathbf{h}\right)\|\right)\|E_{N}\left[\boldsymbol{\xi}_{2}\right]\|$
$\displaystyle+\left(\|\mathbf{Q}\|\|\mathbf{\hat{V}}^{-1}-\mathbf{{V}}^{-1}\|+\|\mathbf{\hat{Q}}-\mathbf{Q}\|\|\mathbf{\hat{V}}^{-1}\|\right)\|E_{N}[\boldsymbol{\xi}_{3}]\|$
$\displaystyle+\|\boldsymbol{\gamma}\|\left\|\mathbf{\hat{h}}-\mathbf{h}\right\|.$

Lemmas 14, 17, 19, 18, 22, 21, and 23 and Theorem 3 bound each term on the
right hand side of the above equation, which leads to

$\boldsymbol{\hat{\theta}}-\boldsymbol{\theta}-(E_{N}-E)[\boldsymbol{\zeta}]=O\left(h_{N}^{L}\pi_{h}\right)+O_{p}\left(\sqrt{\frac{\pi_{h}}{N}}+1\\{\alpha\leq
2\\}\sqrt{\frac{\pi_{h}^{2}}{N^{2}h_{N}^{4}\pi_{h\setminus
0}}}+1\\{\alpha>2\\}\frac{1}{Nh_{N}}\right),$

where we can take a sufficiently small $\epsilon>0$. By (24), we have (26).

Third, we are going to show

$E[\|Var(\boldsymbol{\zeta})^{-1/2}\boldsymbol{\zeta}\|^{4}]=o(N).$ (27)

By Lemmas 23 and 24,

$E[\|Var(\boldsymbol{\zeta})^{-1/2}\boldsymbol{\zeta}\|^{4}]=O\left(1\\{\alpha\leq
2\\}\frac{1}{\pi_{h\setminus 0}^{1+\epsilon}}+1\\{2<\alpha\leq
4\\}\frac{\pi_{h\setminus 0}^{1-\epsilon}}{h_{N}^{4}}+1\\{\alpha>4\\}\right),$

where we can take a sufficiently small $\epsilon>0$. By (24), we have (27).

Now we are going to show the statement of this theorem. By Lemma 24, (25) and
(26),

$\sqrt{N}E_{N}\left[\boldsymbol{\hat{\zeta}}\boldsymbol{\hat{\zeta}}^{\prime}\right]^{-1/2}(\boldsymbol{\hat{\theta}}-\boldsymbol{\theta})=\sqrt{N}Var(\boldsymbol{\zeta})^{-1/2}(E_{N}-E)[\boldsymbol{\zeta}]+o_{p}(1).$

By (27), Assumption 7 (i), and Lyapunov’s central limit theorem,
$\sqrt{N}E_{N}\left[\boldsymbol{\hat{\zeta}}\boldsymbol{\hat{\zeta}}^{\prime}\right]^{-1/2}(\boldsymbol{\hat{\theta}}-\boldsymbol{\theta})\rightarrow_{d}\mathcal{N}\left(0,\mathbf{I}_{p}\right).$
∎

##  Appendix D Auxiliary Results for Theorem 4

We introduce and prove auxiliary lemmas, which we use for the proof of Theorem
4. In all the lemma in this section, we impose the assumptions in the
statement of Theorem 4 as well as $T>p$.

######  Lemma 14.

$\boldsymbol{\Psi}_{1:L}^{-1}\mathbf{h}=O(\frac{\pi_{h}}{h_{N}})$,
$\left\|\frac{h_{N}}{\pi_{h}}\boldsymbol{\Psi}_{1:L}^{-1}\mathbf{h}\right\|$
is bounded away from zero, and
$\boldsymbol{\Psi}_{1:L}^{-1}(\mathbf{\hat{h}}-\mathbf{h})=O_{p}\left(\sqrt{\frac{\pi_{h}}{Nh_{N}^{2}}}\right)$.

###### Proof.

While the definition of $D$ is different, the proof is the same as that of
Lemma 1. ∎

######  Lemma 15.

(i) For every number $l\geq 1$,

$\frac{1}{h_{N}^{l}\pi_{h\setminus 0}}E[D^{l}\boldsymbol{1}\\{D\leq
h_{N}\\}]=\frac{\alpha}{l+\alpha}+o(1)$

as $h_{N}\rightarrow 0$. (ii) For every pair of positive numbers
$l_{1},l_{2}$,
$E\left[D^{l_{1}}\left\|\mathbf{D}_{1:L}^{\prime}\boldsymbol{\Psi}_{1:L}^{-1}\right\|^{l_{2}}\right]=O(\pi_{h\setminus
0}h_{N}^{l_{1}})$. (iii)
$E\left[\|Z\|\|\mathbf{D}_{1:L}^{\prime}\boldsymbol{\Psi}_{1:L}^{-1}\|^{l}\right]=O\left(\pi_{h\setminus
0}\right)$ for every random variable $Z$ in Assumption 7 (iii) and for every
positive integer $l$.

###### Proof.

First, we are going to show the statement in (i). Under Assumption 5,
Karamata’s theorem (Resnick, 2007, Theorem 2.1 and Exercise 2.5) implies

$\frac{1}{h_{N}^{l}\pi_{h\setminus 0}}E\left[\boldsymbol{1}\\{D\leq
h_{N}\\}D^{l}\right]=\frac{E\left[\boldsymbol{1}\\{D^{-1}\geq
h_{N}^{-1}\\}(D^{-1})^{-l}\mid D\neq 0\right]}{(h_{N}^{-1})^{-l}P(D^{-1}\geq
h_{N}^{-1}\mid D\neq 0)}=\frac{\alpha}{l+\alpha}+o(1).$

Therefore, the statement in (i) holds.

Next, we are going to show the statement in (ii). Note that

$D^{l_{1}}\left\|\mathbf{D}_{1:L}^{\prime}\boldsymbol{\Psi}_{1:L}^{-1}\right\|^{l_{2}}=D^{l_{1}}\boldsymbol{1}\\{0<D\leq
h_{N}\\}\left\|\left(\begin{array}[]{cccc}\frac{D}{h_{N}}&\cdots&\frac{D^{L}}{h_{N}^{L}}\end{array}\right)^{\prime}\right\|^{l_{2}}\leq
L^{l_{2}/2}h_{N}^{l_{1}}\boldsymbol{1}\\{0<D\leq h_{N}\\}.$

Therefore,

$E\left[D^{l_{1}}\left\|\mathbf{D}_{1:L}^{\prime}\boldsymbol{\Psi}_{1:L}^{-1}\right\|^{l_{2}}\right]\leq
L^{l_{2}/2}h_{N}^{l_{1}}E[\boldsymbol{1}\\{0<D\leq h_{N}\\}]=O(\pi_{h\setminus
0}h_{N}^{l_{1}}).$

Third, we are going to show the statement in (iii). There is a positive
constant $c$ such that $\sup_{u\in[-c,c]}E[\|Z\|\mid D=u]<\infty$. By the
statement in (ii),

$E\left[\|Z\|\|\mathbf{D}_{1:L}^{\prime}\boldsymbol{\Psi}_{1:L}^{-1}\|^{l}\right]\leq
E\left[\|\mathbf{D}_{1:L}^{\prime}\boldsymbol{\Psi}_{1:L}^{-1}\|^{l}\right]\sup_{u\in[-h_{N},h_{N}]}E[\|Z\|\mid
D=u]=O\left(\pi_{h\setminus 0}\right).$

∎

######  Lemma 16.

(i) For every positive number $l$,

$E\left[\boldsymbol{1}\\{h_{N}<D\\}D^{-l}\right]=O(1)\mbox{ if }l<\alpha$

and

$\frac{h_{N}^{l}}{\pi_{h\setminus
0}}E\left[\boldsymbol{1}\\{h_{N}<D\\}D^{-l}\right]=\frac{\alpha}{l-\alpha}+o(1)\mbox{
if }l>\alpha.$

(ii)$E\left[\boldsymbol{1}\\{D>h_{N}\\}D^{-l}\|Z\|\right]=O\left(1\\{\alpha\leq
l\\}\frac{\pi_{h\setminus 0}^{1-\epsilon}}{h_{N}^{l}}+1\\{\alpha>l\\}\right)$
for every positive number $\epsilon\in(0,1)$ for every random variable $Z$ in
Assumption 7 (ii) and for every positive integer $l$.

###### Proof.

First, we are going to show the statement in (i). Under Assumption 5,
Karamata’s theorem (Resnick, 2007, Theorem 2.1 and Exercise 2.5) implies that

$\frac{h_{N}^{l}}{\pi_{h\setminus
0}}E\left[\boldsymbol{1}\\{D>h_{N}\\}D^{-l}\right]=\frac{E\left[\boldsymbol{1}\\{D^{-1}\leq
h_{N}^{-1}\\}(D^{-1})^{l}\mid D\neq 0\right]}{(h_{N}^{-1})^{l}P(D^{-1}\geq
h_{N}^{-1}\mid D\neq 0)}=\frac{\alpha}{l-\alpha}+o(1).$

Next, we are going to show the statement in (ii). There is a positive constant
$c$ such that $\sup_{u\in[-c,c]}E[\|Z\|\mid D=u]<\infty$. Note that

$\displaystyle\left\|E\left[\boldsymbol{1}\\{|D|>h_{N}\\}|D|^{-l}\|Z\|\right]\right\|\leq$
$\displaystyle\left\|E\left[\boldsymbol{1}\\{|D|>c\\}|D|^{-l}\|Z\|\right]\right\|+\left\|E\left[\boldsymbol{1}\\{h_{N}<|D|\leq
c\\}|D|^{-l}\|Z\|\right]\right\|$ $\displaystyle\leq$ $\displaystyle
c^{-l}E\left[\|Z\|\right]+E\left[\boldsymbol{1}\\{h_{N}<|D|\\}|D|^{-l}\right]\sup_{u\in[-c,c]}E[\|Z\|\mid
D=u].$

By Hölder’s inequality and the statement in (i),

$E\left[\boldsymbol{1}\\{|D|>h_{N}\\}|D|^{-l}\right]=\begin{cases}O(1)&\mbox{
if }\alpha>l\\\ O\left(\frac{\pi_{h\setminus 0}}{h_{N}^{l}}\right)&\mbox{ if
}\alpha<l\\\
O(E\left[\boldsymbol{1}\\{|D|>h_{N}\\}|D|^{-l/{1-\epsilon}}\right]^{1-\epsilon})=O\left(\frac{\pi_{h\setminus
0}}{h_{N}^{l/{1-\epsilon}}}\right)^{1-\epsilon}=O\left(\frac{\pi_{h\setminus
0}^{1-\epsilon}}{h_{N}^{l}}\right)&\mbox{ if }\alpha=l,\end{cases}$

which implies the statement in (ii). ∎

######  Lemma 17.

The minimum eigenvalue of $\boldsymbol{\Pi}$ is bounded away from zero, the
maximum eigenvalue of $\boldsymbol{\Pi}$ is bounded, and
$\boldsymbol{\hat{\Pi}}^{-1}-\boldsymbol{\Pi}^{-1}=O_{p}\left(\frac{1}{\sqrt{N\pi_{h\setminus
0}}}\right)$.

###### Proof.

While the definition of $D$ is different, the proof is the same as that of
Lemma 4. ∎

######  Lemma 18.

$\mathbf{V}^{-1}=O(1)$, and
$\mathbf{\hat{V}}^{-1}-\mathbf{V}^{-1}=O_{p}\left(\sqrt{\frac{1}{N}}\right)$.

###### Proof.

By Assumption 7 (iv), the minimum eigenvalue of $\mathbf{V}$ is bounded away
from zero, which implies $\mathbf{V}^{-1}=O(1)$. Since $\mathbf{V}^{-1}=O(1)$,
it suffices to show
$\mathbf{\hat{V}}-\mathbf{V}=O_{p}\left(\sqrt{\frac{1}{N}}\right)$ to show
$\mathbf{\hat{V}}^{-1}-\mathbf{V}^{-1}=O_{p}\left(\sqrt{\frac{1}{N}}\right)$.
Let $\mathbf{c}$ be any $(T-1)p$-dimensional column vector with $\|c\|=1$.
Together with Assumption 7 (i), we have

$E\left[\left\|(\mathbf{\hat{V}}-\mathbf{V})\mathbf{c}\right\|^{2}\right]^{1/2}=E\left[\left\|(E_{N}-E)\left[\boldsymbol{1}\\{D>h_{N}\\}\mathbf{W}^{\prime}\mathbf{M}_{\mathbf{X}}\mathbf{W}\mathbf{c}\right]\right\|^{2}\right]^{1/2}\leq\frac{1}{\sqrt{N}}E\left[\|\mathbf{W}\|^{4}\|\mathbf{M}_{\mathbf{X}}\|^{2}\right]^{1/2}.$

Since $\|\mathbf{M}_{\mathbf{X}}\|\leq 1$ and
$E\left[\|\mathbf{W}\|^{4}\right]<\infty$, we have
$\mathbf{\hat{V}}-\mathbf{V}=O_{p}\left(\sqrt{\frac{1}{N}}\right)$. ∎

######  Lemma 19.

For every $\epsilon\in(0,1)$, $\mathbf{Q}=O(1\\{\alpha\leq
1\\}\frac{\pi_{h}^{1-\epsilon}}{h_{N}}+1\\{\alpha>1\\})$, and
$\mathbf{\hat{Q}}-\mathbf{Q}=O_{p}\Big{(}1\\{\alpha\leq
2\\}\sqrt{\frac{\pi_{h}^{2-\epsilon}}{Nh_{N}^{2}\pi_{h\setminus
0}}}+1\\{\alpha>2\\}\sqrt{\frac{1}{N}}\Big{)}$.

###### Proof.

Let $\mathbf{c}$ be any $(T-1)p$-dimensional column vector with
$\|\mathbf{c}\|=1$. The the first statement follows from Lemmas 14, 15 and 16
and

$\displaystyle\mathbf{Q}\mathbf{c}$ $\displaystyle=$
$\displaystyle\mathbf{R}\mathbf{c}+E\left[\boldsymbol{1}\\{D>h_{N}\\}D^{-1}(\mathbf{X}^{\prime}\mathbf{X})^{\ast}\mathbf{X}^{\prime}\mathbf{W}\mathbf{c}\right]+\frac{1}{\pi_{h\setminus
0}}E\left[(\mathbf{X}^{\prime}\mathbf{X})^{\ast}\mathbf{X}^{\prime}\mathbf{W}\mathbf{c}\mathbf{D}_{1:L}^{\prime}\boldsymbol{\Psi}_{1:L}^{-1}\right]\boldsymbol{\Pi}^{-1}\boldsymbol{\Psi}_{1:L}^{-1}\mathbf{h}$
$\displaystyle=$ $\displaystyle O\left(1\\{\alpha\leq
1\\}\frac{\pi_{h}^{1-\epsilon}}{h_{N}}+1\\{\alpha>1\\}\right).$

Assumption 7 (i) and Lemmas 15 and 16 imply

$\|E\left[(\mathbf{X}^{\prime}\mathbf{X})^{\ast}\mathbf{X}^{\prime}\mathbf{W}\mathbf{c}\mathbf{D}_{1:L}^{\prime}\boldsymbol{\Psi}_{1:L}^{-1}\right]\|\leq
E\left[\|(\mathbf{X}^{\prime}\mathbf{X})^{\ast}\mathbf{X}^{\prime}\mathbf{W}\mathbf{c}\|\|\mathbf{D}_{1:L}^{\prime}\boldsymbol{\Psi}_{1:L}^{-1}\|\right]=O\left(\pi_{h\setminus
0}\right),$
$E\left[\left\|(E_{N}-E)\left[(\mathbf{X}^{\prime}\mathbf{X})^{\ast}\mathbf{X}^{\prime}\mathbf{W}\mathbf{c}\mathbf{D}_{1:L}^{\prime}\boldsymbol{\Psi}_{1:L}^{-1}\right]\right\|^{2}\right]^{1/2}\leq\frac{1}{\sqrt{N}}E\left[\left\|(\mathbf{X}^{\prime}\mathbf{X})^{\ast}\mathbf{X}^{\prime}\mathbf{W}\mathbf{c}\right\|^{2}\left\|\mathbf{D}_{1:L}^{\prime}\boldsymbol{\Psi}_{1:L}^{-1}\right\|^{2}\right]^{1/2}=O\left(\sqrt{\frac{\pi_{h\setminus
0}}{N}}\right),$

and

$\displaystyle
E\left[\left\|(E_{N}-E)\left[(\mathbf{X}^{\prime}\mathbf{X})^{\ast}\mathbf{X}^{\prime}\mathbf{W}\boldsymbol{1}\\{D>h_{N}\\}D^{-1}\right]\right\|^{2}\right]^{1/2}\leq$
$\displaystyle\frac{1}{\sqrt{N}}E\left[\left\|(\mathbf{X}^{\prime}\mathbf{X})^{\ast}\mathbf{X}^{\prime}\mathbf{W}\right\|^{2}\boldsymbol{1}\\{D>h_{N}\\}D^{-2}\right]^{1/2}$
$\displaystyle=$ $\displaystyle O\left(1\\{\alpha\leq
2\\}\sqrt{\frac{\pi_{h\setminus
0}^{1-\epsilon}}{Nh_{N}^{2}}}+1\\{\alpha>2\\}\sqrt{\frac{1}{N}}\right).$

Since

$\displaystyle(\mathbf{\hat{Q}}-\mathbf{Q})\mathbf{c}=$
$\displaystyle\frac{1}{\pi_{h\setminus
0}}E\left[(\mathbf{X}^{\prime}\mathbf{X})^{\ast}\mathbf{X}^{\prime}\mathbf{W}\mathbf{c}\mathbf{D}_{1:L}^{\prime}\boldsymbol{\Psi}_{1:L}^{-1}\right](\boldsymbol{\hat{\Pi}}^{-1}-\boldsymbol{\Pi}^{-1})\boldsymbol{\Psi}_{1:L}^{-1}\mathbf{h}$
$\displaystyle+\frac{1}{\pi_{h\setminus
0}}(E_{N}-E)\left[(\mathbf{X}^{\prime}\mathbf{X})^{\ast}\mathbf{X}^{\prime}\mathbf{W}\mathbf{c}\mathbf{D}_{1:L}^{\prime}\boldsymbol{\Psi}_{1:L}^{-1}\right]\boldsymbol{\hat{\Pi}}^{-1}\boldsymbol{\Psi}_{1:L}^{-1}\mathbf{h}$
$\displaystyle+(E_{N}-E)\left[(\mathbf{X}^{\prime}\mathbf{X})^{\ast}\mathbf{X}^{\prime}\mathbf{W}\boldsymbol{1}\\{D>h_{N}\\}D^{-1}\right]\mathbf{c},$

together with Lemmas 14 and 17, we have

$(\mathbf{\hat{Q}}-\mathbf{Q})\mathbf{c}=O_{p}\left(1\\{\alpha\leq
2\\}\sqrt{\frac{\pi_{h}^{2-\epsilon}}{Nh_{N}^{2}\pi_{h\setminus
0}}}+1\\{\alpha>2\\}\sqrt{\frac{1}{N}}\right).$

∎

######  Lemma 20.

For every $\epsilon\in(0,1)$, the following statements hold: (i)
$E[\boldsymbol{\xi}_{1}]=0$. (ii)
$E[\|\boldsymbol{\xi}_{1}\|^{2}]=O(1\\{\alpha\leq 2\\}\frac{\pi_{h\setminus
0}^{1-\epsilon}}{h_{N}^{2}}+1\\{\alpha>2\\})$ and
$E[\|\boldsymbol{\xi}_{1}\|^{4}]=O(1\\{\alpha\leq 4\\}\frac{\pi_{h\setminus
0}^{1-\epsilon}}{h_{N}^{4}}+1\\{\alpha>4\\})$. (iii)
$(E_{N}-E)[\boldsymbol{\xi}_{1}]=O_{p}(1\\{\alpha\leq
2\\}\sqrt{\frac{\pi_{h\setminus
0}^{1-\epsilon}}{Nh_{N}^{2}}}+1\\{\alpha>2\\}\sqrt{\frac{1}{N}})$.

###### Proof.

The statement in (i) follows from the definition of $\boldsymbol{\xi}_{1}$.
The statement in (ii) follows from Lemma 16. The statement in (iii) follows
from Assumption 7 (i) and the statement in (ii). ∎

######  Lemma 21.

(i) $E[\boldsymbol{\xi}_{2}]=O(h_{N}^{L+1}\pi_{h\setminus 0})$. (ii)
$E[\|\boldsymbol{\xi}_{2}\|^{2}]=O\left(\pi_{h\setminus 0}\right)$ and
$E[\|\boldsymbol{\xi}_{2}\|^{4}]=O\left(\pi_{h\setminus 0}\right)$. (iii)
$(E_{N}-E)[\boldsymbol{\xi}_{2}]=O(\sqrt{\frac{\pi_{h\setminus 0}}{N}})$.

###### Proof.

By Assumption 6 (i) and $\mathbf{m}(0)=0$, there is a constant ${C}$ such that
$\left\|\mathbf{m}(u)-\sum_{l=1}^{L}\mathbf{m}^{(l)}(0)\frac{u^{l}}{l!}\right\|\leq{C}|u|^{L+1}$
for every $u\in[-h_{N},h_{N}]$ when $h_{N}$ is sufficiently small. Therefore,
for every $u\in[-h_{N},h_{N}]$, we have

$\displaystyle\left\|E\left[\boldsymbol{\xi}_{2}\right]\right\|=$
$\displaystyle\left\|E\left[\left(\mathbf{m}(D)-\sum_{l=1}^{L}\mathbf{m}^{(l)}(0)\frac{D^{l}}{l!}\right)\mathbf{D}_{1:L}^{\prime}\boldsymbol{\Psi}_{1:L}^{-1}\right]\right\|$
$\displaystyle\leq$ $\displaystyle
E\left[\left\|\mathbf{m}(D)-\sum_{l=1}^{L}\mathbf{m}^{(l)}(0)\frac{D^{l}}{l!}\right\|\left\|\mathbf{D}_{1:L}^{\prime}\boldsymbol{\Psi}_{1:L}^{-1}\right\|\right]$
$\displaystyle\leq$
$\displaystyle{C}E\left[D^{L+1}\left\|\mathbf{D}_{1:L}^{\prime}\boldsymbol{\Psi}_{1:L}^{-1}\right\|\right]$
$\displaystyle=$ $\displaystyle O(\pi_{h\setminus 0}h_{N}^{L+1}),$

where the last equality follows from Lemma 15. This shows the statement in
(i). The statement in (ii) follows from Lemma 15. The statement in (iii)
follows from Assumption 7 (i) and the statement in (ii). ∎

######  Lemma 22.

For every $\epsilon\in(0,1)$, the following statements hold: (i)
$E\left[\boldsymbol{\xi}_{3}\right]=0$. (ii)
$E[\|\boldsymbol{\xi}_{3}\|^{2}]=O\left(1\right)$ and
$E[\|\boldsymbol{\xi}_{3}\|^{4}]=O\left(1\right)$. (iii)
$(E_{N}-E)\left[\boldsymbol{\xi}_{3}\right]=O_{p}\left(\sqrt{\frac{1}{N}}\right)$.

###### Proof.

The first statement follows from (11). Since
$\|\boldsymbol{\xi}_{3}\|\leq\|\mathbf{W}\|\|\mathbf{M}_{\mathbf{X}}\|\|\mathbf{Y}-\mathbf{W}\boldsymbol{\delta}\|$
and $\|\mathbf{M}_{\mathbf{X}}\|\leq 1$, we have

$E[\|\boldsymbol{\xi}_{3}\|^{l}]^{1/l}\leq
E[\|\mathbf{W}\|^{l}\|\mathbf{Y}-\mathbf{W}\boldsymbol{\delta}\|^{l}]^{1/l}$

for $l=2,4$, which implies the statement in (ii). The statement in (iii)
follows from Assumption 7 (i) and the statement in (i). ∎

######  Lemma 23.

For every $\epsilon\in(0,1)$, the following statements hold: (i)
$E[\boldsymbol{\zeta}]=O(h_{N}^{L}\pi_{h})$. (ii)
$E[\|\boldsymbol{\zeta}\|^{2}]=O(1\\{\alpha\leq
1\\}\frac{\pi_{h}^{2-\epsilon}}{h_{N}^{2}\pi_{h\setminus 0}}+1\\{1<\alpha\leq
2\\}\frac{\pi_{h\setminus 0}^{1-\epsilon}}{h_{N}^{2}}+1\\{\alpha>2\\})$ and
$E[\|\boldsymbol{\zeta}\|^{4}]=O(1\\{\alpha\leq
1\\}\frac{\pi_{h}^{4-\epsilon}}{h_{N}^{4}\pi_{h\setminus
0}^{3}}+1\\{1<\alpha\leq 4\\}\frac{\pi_{h\setminus
0}^{1-\epsilon}}{h_{N}^{4}}+1\\{\alpha>4\\})$. (iii)
$(E_{N}-E)[\boldsymbol{\zeta}\boldsymbol{\zeta}^{\prime}]=O_{p}\left(1\\{\alpha\leq
1\\}\sqrt{\frac{\pi_{h}^{4-\epsilon}}{Nh_{N}^{4}\pi_{h\setminus
0}^{3}}}+1\\{1<\alpha\leq 4\\}\sqrt{\frac{\pi_{h\setminus
0}^{1-\epsilon}}{Nh_{N}^{4}}}+1\\{\alpha>4\\}\sqrt{\frac{1}{N}}\right)$.

###### Proof.

By Lemmas 14, 17, and 21, we have

$E[\boldsymbol{\zeta}]=\frac{1}{\pi_{h\setminus
0}}E[\boldsymbol{\xi}_{2}]\boldsymbol{\Pi}^{-1}\boldsymbol{\Psi}_{1:L}^{-1}\mathbf{h}=O(h_{N}^{L}\pi_{h}),$

which proves the statement in (i). Note that

$E[\|\boldsymbol{\zeta}\|^{l}]^{1/l}\leq
E[\|\boldsymbol{\xi}_{1}\|^{l}]^{1/l}+\frac{1}{\pi_{h\setminus
0}}E[\|\boldsymbol{\xi}_{2}\|^{l}]^{1/l}\|\boldsymbol{\Pi}^{-1}\|\|\boldsymbol{\Psi}_{1:L}^{-1}\mathbf{h}\|+\|\mathbf{Q}\|\|\mathbf{{V}}^{-1}\|E[\|\boldsymbol{\xi}_{3}\|^{l}]^{1/l}$

for $l=2,4$. By Lemmas 19, 20, 21, and 22, we have

$E[\|\boldsymbol{\zeta}\|^{l}]=O\left(1\\{\alpha\leq
1\\}\frac{\pi_{h}^{l-\epsilon}}{h_{N}^{l}\pi_{h\setminus
0}^{l-1}}+1\\{1<\alpha\leq l\\}\frac{\pi_{h\setminus
0}^{1-\epsilon}}{h_{N}^{l}}+1\\{\alpha>l\\}\right),$

which leads to the statement in (ii). The statement in (iii) follows from
Assumption 7 (i) and the statement in (ii). ∎

######  Lemma 24.

The minimum eigenvalue of $\left(1\\{\alpha\leq
2\\}\frac{h_{N}^{2}\pi_{h\setminus
0}}{\pi_{h}^{2}}+1\\{\alpha>2\\}\right)Var(\boldsymbol{\zeta})$ is bounded
away from zero.

###### Proof.

Let $\mathbf{c}$ be any $p$-dimensional column vector with $\|\mathbf{c}\|=1$.
In $\boldsymbol{\zeta}$, the terms
$\left(\boldsymbol{\xi}_{1}+\mathbf{Q}\mathbf{{V}}^{-1}\boldsymbol{\xi}_{3}\right)$
and $\frac{1}{\pi_{h\setminus
0}}\boldsymbol{\xi}_{2}\boldsymbol{\Pi}^{-1}\boldsymbol{\Psi}_{1:L}^{-1}\mathbf{h}$
are independent. Therefore,

$Var(\mathbf{c}^{\prime}\boldsymbol{\zeta})\geq\max\left\\{Var(\mathbf{c}^{\prime}\left(\boldsymbol{\xi}_{1}+\mathbf{Q}\mathbf{{V}}^{-1}\boldsymbol{\xi}_{3}\right)),\frac{1}{\pi_{h\setminus
0}^{2}}Var(\mathbf{c}^{\prime}\boldsymbol{\xi}_{2}\boldsymbol{\Pi}^{-1}\boldsymbol{\Psi}_{1:L}^{-1}\mathbf{h})\right\\}.$

First, we are going to show that $\frac{h_{N}^{2}}{\pi_{h}^{2}\pi_{h\setminus
0}}Var(\mathbf{c}^{\prime}\boldsymbol{\xi}_{2}\boldsymbol{\Pi}^{-1}\boldsymbol{\Psi}_{1:L}^{-1}\mathbf{h})$
is bounded away from zero when $\alpha\leq 2$. By Assumptions 6 (ii) and 7
(iv), there is a positive constant $C$ such that
$\mathbf{c}^{\prime}Var\left((\mathbf{X}^{\prime}\mathbf{X})^{\ast}\mathbf{X}^{\prime}(\mathbf{Y}-\mathbf{W}\boldsymbol{\delta})\mid
D=u\right)\mathbf{c}\geq C$ for every $u\in[-h_{N},h_{N}]$ for sufficiently
small $h_{N}$. By the conditional variance identity,

$Var(\mathbf{c}^{\prime}\boldsymbol{\xi}_{2}\boldsymbol{\Pi}^{-1}\boldsymbol{\Psi}_{1:L}^{-1}\mathbf{h})\geq
E[Var((\boldsymbol{\xi}_{2}\boldsymbol{\Pi}^{-1}\boldsymbol{\Psi}_{1:L}^{-1}\mathbf{h})^{\prime}\mathbf{c}\mid
D)].$

Since

$\displaystyle
Var((\boldsymbol{\xi}_{2}\boldsymbol{\Pi}^{-1}\boldsymbol{\Psi}_{1:L}^{-1}\mathbf{h})^{\prime}\mathbf{c}\mid
D)=$
$\displaystyle\left(\boldsymbol{\Psi}_{1:L}^{-1}\mathbf{h}\right)^{\prime}\boldsymbol{\Pi}^{-1}(\mathbf{c}^{\prime}Var((\mathbf{X}^{\prime}\mathbf{X})^{\ast}\mathbf{X}^{\prime}(\mathbf{Y}-\mathbf{W}\boldsymbol{\delta})\mid
D)\mathbf{c})\boldsymbol{\Psi}_{1:L}^{-1}\mathbf{D}_{1:L}\mathbf{D}_{1:L}^{\prime}\boldsymbol{\Psi}_{1:L}^{-1}\boldsymbol{\Pi}^{-1}\boldsymbol{\Psi}_{1:L}^{-1}\mathbf{h}$
$\displaystyle\geq$ $\displaystyle
C\left(\boldsymbol{\Psi}_{1:L}^{-1}\mathbf{h}\right)^{\prime}\boldsymbol{\Pi}^{-1}\boldsymbol{\Psi}_{1:L}^{-1}\mathbf{D}_{1:L}\mathbf{D}_{1:L}^{\prime}\boldsymbol{\Psi}_{1:L}^{-1}\boldsymbol{\Pi}^{-1}\boldsymbol{\Psi}_{1:L}^{-1}\mathbf{h},$

we have

$Var(\mathbf{c}^{\prime}\boldsymbol{\xi}_{2}\boldsymbol{\Pi}^{-1}\boldsymbol{\Psi}_{1:L}^{-1}\mathbf{h})\geq\pi_{h\setminus
0}C\left(\boldsymbol{\Psi}_{1:L}^{-1}\mathbf{h}\right)^{\prime}\boldsymbol{\Pi}^{-1}\boldsymbol{\Psi}_{1:L}^{-1}\mathbf{h}.$

By Lemmas 14 and 17, $\frac{h_{N}^{2}}{\pi_{h}^{2}\pi_{h\setminus
0}}Var(\mathbf{c}^{\prime}\boldsymbol{\xi}_{2}\boldsymbol{\Pi}^{-1}\boldsymbol{\Psi}_{1:L}^{-1}\mathbf{h})$
is bounded away from zero.

Second, we are going to show that
$Var(\mathbf{c}^{\prime}(\boldsymbol{\xi}_{1}+\mathbf{Q}\mathbf{{V}}^{-1}\boldsymbol{\xi}_{3}))$
is bounded away from zero when $\alpha>2$. By Assumption 7 (iv),
$Var(\mathbf{c}^{\prime}\boldsymbol{\eta})$ is bounded away from zero. Note
that

$\displaystyle
E\left[\left\|\mathbf{c}^{\prime}\left(\boldsymbol{\eta}-(\boldsymbol{\xi}_{1}+\mathbf{Q}\mathbf{{V}}^{-1}\boldsymbol{\xi}_{3})\right)\right\|^{2}\right]^{1/2}=$
$\displaystyle E\left[\left\|\frac{1}{\pi_{h\setminus
0}}\boldsymbol{\xi}_{3}^{\prime}\mathbf{{V}}^{-1}E\left[\mathbf{W}^{\prime}\mathbf{X}(\mathbf{X}^{\prime}\mathbf{X})^{\ast}\mathbf{c}\mathbf{D}_{1:L}^{\prime}\boldsymbol{\Psi}_{1:L}^{-1}\right]\boldsymbol{\Pi}^{-1}\boldsymbol{\Psi}_{1:L}^{-1}\mathbf{h}\right\|^{2}\right]^{1/2}$
$\displaystyle\leq$ $\displaystyle\frac{1}{\pi_{h\setminus
0}}E\left[\left\|\boldsymbol{\xi}_{3}\right\|^{2}\right]^{1/2}\|\mathbf{{V}}^{-1}\|E\left[\|\mathbf{W}^{\prime}\mathbf{X}(\mathbf{X}^{\prime}\mathbf{X})^{\ast}\|\|\mathbf{D}_{1:L}^{\prime}\boldsymbol{\Psi}_{1:L}^{-1}\|\right]\|\boldsymbol{\Pi}^{-1}\|\|\boldsymbol{\Psi}_{1:L}^{-1}\mathbf{h}\|$
$\displaystyle=$ $\displaystyle O\left(\frac{\pi_{h}}{h_{N}}\right)$
$\displaystyle=$ $\displaystyle o(1),$

where the last equality uses $\alpha>1$. Therefore,
$Var(\mathbf{c}^{\prime}\left(\boldsymbol{\xi}_{1}+\mathbf{Q}\mathbf{{V}}^{-1}\boldsymbol{\xi}_{3}\right))$
is bounded away from zero. ∎

######  Lemma 25.

$E_{N}[\|\boldsymbol{\hat{\zeta}}-\boldsymbol{\zeta}\|^{2}]^{1/2}=o(1)+O_{p}\left(\sqrt{\frac{\pi_{h}^{2}}{Nh_{N}^{2}\pi_{h\setminus
0}^{2}}}\right)$.

###### Proof.

By Lemmas 14, 15, 16, 17, 18, 19, 21, 22, 23, and 25, we have

$\displaystyle\|\boldsymbol{\hat{\delta}}-\boldsymbol{\delta}\|\leq$
$\displaystyle\|\mathbf{\hat{V}}^{-1}\|\|E[\boldsymbol{\xi}_{3}]\|+\|\mathbf{\hat{V}}^{-1}\|\|(E_{N}-E)[\boldsymbol{\xi}_{3}]\|=O_{p}\left(\sqrt{\frac{1}{N}}\right)$
$\displaystyle\|(\boldsymbol{\hat{\gamma}}-\boldsymbol{\gamma})\boldsymbol{\Psi}_{1:L}\|\leq$
$\displaystyle\frac{1}{\pi_{h\setminus
0}}\|E\left[\boldsymbol{\xi}_{2}\right]\|\|\boldsymbol{\hat{\Pi}}^{-1}\|+\frac{1}{\pi_{h\setminus
0}}\|(E_{N}-E)\left[\boldsymbol{\xi}_{2}\right]\|\|\boldsymbol{\hat{\Pi}}^{-1}\|$
$\displaystyle+\frac{1}{\pi_{h\setminus
0}}\|\boldsymbol{\hat{\delta}}-\boldsymbol{\delta}\|E_{N}\left[\|(\mathbf{X}^{\prime}\mathbf{X})^{\ast}\mathbf{X}^{\prime}\mathbf{W}\|\|\mathbf{D}_{1:L}^{\prime}\boldsymbol{\Psi}_{1:L}^{-1}\|\right]\|\boldsymbol{\hat{\Pi}}^{-1}\|$
$\displaystyle=$ $\displaystyle
O(h_{N}^{L+1})+O_{p}\left(\sqrt{\frac{1}{N\pi_{h\setminus 0}}}\right)$
$\displaystyle
E_{N}[\|\boldsymbol{\hat{\xi}}_{1}-\boldsymbol{\xi}_{1}\|^{2}]^{1/2}\leq$
$\displaystyle\|E_{N}[\boldsymbol{\xi}_{1}]\|+E_{N}[\boldsymbol{1}\\{D>h_{N}\\}D^{-2}\|(\mathbf{X}^{\prime}\mathbf{X})^{\ast}\mathbf{X}^{\prime}\mathbf{W}\|^{2}]^{1/2}\|\boldsymbol{\hat{\delta}}-\boldsymbol{\delta}\|=O_{p}\left(\sqrt{\frac{\pi_{h}^{2}}{Nh_{N}^{2}\pi_{h\setminus
0}^{2}}}\right),$ $\displaystyle
E_{N}[\|\boldsymbol{\hat{\xi}}_{2}-\boldsymbol{\xi}_{2}\|^{2}]^{1/2}\leq$
$\displaystyle
E_{N}[\|(\mathbf{X}^{\prime}\mathbf{X})^{\ast}\mathbf{X}^{\prime}\mathbf{W}\|^{2}\|\mathbf{D}_{1:L}^{\prime}\boldsymbol{\Psi}_{1:L}^{-1}\|^{2}]^{1/2}\|\boldsymbol{\hat{\delta}}-\boldsymbol{\delta}\|+\|(\boldsymbol{\hat{\gamma}}-\boldsymbol{\gamma})\boldsymbol{\Psi}_{1:L}\|E_{N}[\|\mathbf{D}_{1:L}^{\prime}\boldsymbol{\Psi}_{1:L}^{-1}\|^{4}]^{1/2}$
$\displaystyle=$ $\displaystyle
o(h_{N}^{L+1})+O_{p}\left(\frac{1}{\sqrt{N}}\right),$ $\displaystyle
E_{N}[\|\boldsymbol{\hat{\xi}}_{3}-\boldsymbol{\xi}_{3}\|^{2}]^{1/2}\leq$
$\displaystyle
E_{N}[\|\mathbf{W}\|^{4}\|\mathbf{M}_{\mathbf{X}}\|^{2}]^{1/2}\|\boldsymbol{\hat{\delta}}-\boldsymbol{\delta}\|=O_{p}\left(\frac{1}{\sqrt{N}}\right),\mbox{
and }$ $\displaystyle
E_{N}[\|\boldsymbol{\hat{\zeta}}-\boldsymbol{\zeta}\|^{2}]^{1/2}\leq$
$\displaystyle
E_{N}[\|\boldsymbol{\hat{\xi}}_{1}-\boldsymbol{\xi}_{1}\|^{2}]^{1/2}+\frac{1}{\pi_{h\setminus
0}}E_{N}[\|(\boldsymbol{\hat{\xi}}_{2}-\boldsymbol{\xi}_{2})\|^{2}]^{1/2}\|\boldsymbol{\hat{\Pi}}^{-1}\|\|\boldsymbol{\Psi}_{1:L}^{-1}\mathbf{\hat{h}}\|$
$\displaystyle+\left\|\mathbf{\hat{Q}}\right\|\|\mathbf{\hat{V}}^{-1}\|E_{N}[\|\boldsymbol{\hat{\xi}}_{3}-\boldsymbol{\xi}_{3}\|^{2}]^{1/2}$
$\displaystyle+\frac{1}{\pi_{h\setminus
0}}\left(\|\boldsymbol{\hat{\Pi}}^{-1}-\boldsymbol{\Pi}^{-1}\|\|\boldsymbol{\Psi}_{1:L}^{-1}\mathbf{\hat{h}}\|+\|\boldsymbol{\Pi}^{-1}\|\|\boldsymbol{\Psi}_{1:L}^{-1}(\mathbf{\hat{h}}-\mathbf{h})\|\right)E_{N}[\|\boldsymbol{\xi}_{2}\|^{2}]^{1/2}$
$\displaystyle+\left(\left\|\mathbf{\hat{Q}}\right\|\|\mathbf{\hat{V}}^{-1}-\mathbf{{V}}^{-1}\|+\left\|\mathbf{\hat{Q}}-\mathbf{Q}\right\|\|\mathbf{{V}}^{-1}\|\right)E_{N}[\|\boldsymbol{\xi}_{3}\|^{2}]^{1/2}$
$\displaystyle=$ $\displaystyle
o(1)+O_{p}\left(\sqrt{\frac{\pi_{h}^{2}}{Nh_{N}^{2}\pi_{h\setminus
0}^{2}}}\right).$

Therefore, the statement of the lemma follows. ∎


Fast Solver for Diffusive Transport Times on Dynamic NetworksL. Elam, M. C.
Quiñones-Frías, Y. Zhang, A. A. Rodal, and T. G. Fai

# Fast Solver for Diffusive Transport Times on Dynamic Intracellular
Networks††thanks:  We acknowledge funding from the Brandeis NSF MRSEC
Bioinspired Soft Materials (NSF-DMR2011846) and NIH R01 NS116375 (AAR and
TGF), and NSF grants DMS-1913093 and MCB-2213583 (TGF).

Lachlan Elam Department of Mathematics, Brandeis University, Waltham, MA .
Mónica C. Quiñones-Frías Department of Biology, Brandeis University, Waltham,
MA. Ying Zhang22footnotemark: 2 Avital A. Rodal33footnotemark: 3 Thomas G.
Fai22footnotemark: 2<EMAIL_ADDRESS>
###### Abstract

The transport of particles in cells is influenced by the properties of
intracellular networks they traverse while searching for localized target
regions or reaction partners. Moreover, given the rapid turnover in many
intracellular structures, it is crucial to understand how temporal changes in
the network structure affect diffusive transport. In this work, we use network
theory to characterize complex intracellular biological environments across
scales. We develop an efficient computational method to compute the mean first
passage times for simulating a particle diffusing along two-dimensional planar
networks extracted from fluorescence microscopy imaging. We first benchmark
this methodology in the context of synthetic networks, and subsequently apply
it to live-cell data from endoplasmic reticulum tubular networks.

###### keywords:

intracellular transport, mean-first passage time, Sherman-Morrison formula,
dynamic networks, endoplasmic reticulum

##  1 Introduction

Previous experimental and theoretical studies reveal how the diffusive
exploration of objects relates to the morphological properties of
intracellular environments [19, 9, 2, 7]. One particularly important category
of intracellular diffusive transport emerges on network-like structures used
to describe the morphologies of subcellular organelles [17]. Lipids and
proteins are synthesized in the ER, and their rate and direction of transport
from their site of synthesis to other points in the ER network, where they can
be delivered to other organelles, remain poorly understood. Objects traversing
along the edges of a network exhibit a broad range of transport behaviors
dependent on the properties of the network. For example, previous work has
demonstrated the importance of the edge-length distribution in overall search
times on networks through analyzing the mean first passage time (MFPT) on
planar networks [4]. Motivated by the rapid turnover observed in subcellular
structures, in this work we develop an efficient computational method to
compute how changes in network morphology affect the diffusive MFPT.

Changes in the structure of confined networks have been shown to have a
significant impact on transport in a broad variety of biological problems,
including cytoskeletal intracellular transport [1] and cargo trafficking along
the tubules of the endoplasmic reticulum (ER) [11]. In this work, we use a
coarse-grained model of random walks on networks to simulate the diffusion on
networks. By developing an efficient approach based on the Sherman-Morrison
formula, we are able to calculate how changes in the network affect transport.
In particular, we focus on ER networks obtained from experiments and apply the
method to show how the properties of the network influence the MFPT of
particles diffusing along ER tubules. We find that changes in the ER network
structure within dense and highly-connected regions (e.g. in a neighborhood
around the nucleus) tend to speed up the diffusive search of a particle. On
the other hand, the loss of tubules that bridge together different components
of the network tends to slow down the search process.

An open-source repository with implementations of the methods described here
is available at https://github.com/Lanlach1/SIAM-Diffusive-Transport-Time-
Solver.

##  2 Mathematical Model

###  2.1 Transition probabilities along edges of a static network

The networks that we consider are a system of nodes connected by edges with
each node having $N_{i}$ ($N_{i}\geq 1$) connections and particles searching
along a network via transiting from nodes. Assuming particles undergo an
unbiased random walk along the edges of the network, the transition
probabilities satisfy [8, 4, 12]:

(1) $P_{i,j}=\frac{1}{l_{j}}\frac{1}{\sum_{j=1}^{N_{i}}\frac{1}{l_{j}}},$

where $P_{i,j}$, $j=1,2,\dots,N_{i}$ denote the transition probabilities from
$X_{i}$ to its $j^{\text{th}}$ connected node, $X_{j}$, and $l_{j}$ is the
length of the edge between $X_{i}$ and $X_{j}$.

Similarly, we define $Q_{i}$ as the average time spent at node $i$ prior to
stepping to one of its neighbors, which from [4] satisfies

(2)
$Q_{i}=\frac{1}{2D}\frac{\sum_{j=1}^{N_{i}}l_{j}}{\sum_{j=1}^{N_{i}}l_{j}^{-1}},$

where $D$ is the diffusion coefficient. In this work, we take $D=1$ µm2/s
throughout, which gives an order of magnitude estimate for typical proteins
diffusing in cells [10].

We define the stationary distribution $\mathbf{z}\in\mathbb{R}^{n}$ from the
transition matrix $P$ such that $\mathbf{z}^{T}P=\mathbf{z}^{T}$. As defined,
$P$ is a row stochastic matrix, therefore we may compute $\mathbf{z}$ by
solving for the left eigenvector of $P$ associated with eigenvalue 1. Each
node $i$ will have steady state probability $z_{i}$.

###  2.2 Mean first passage time

We quantify the number of nodes visited on average during diffusion between
any two nodes across the network using the mean first passage time (MFPT),
which can be computed from the transition matrix $P$ with entries given by Eq.
1. It follows from the Markov property of the random walk that the transition
probability of reaching node $X_{i}$ from $X_{j}$ in two steps is given by

$\displaystyle
P_{i,j}^{2}=\sum_{k}P_{ik}P_{kj}=\mathbbm{P}\left[\text{Transition to
}X_{i}|\text{Current Node}=X_{j}\,\text{and}\,2\,\text{steps}\right].$

We may further generalize this property to obtain the transition probability
of reaching node $i$ from node $j$ in 3 steps:

$\displaystyle
P_{i,j}^{3}=\sum_{h}P_{ih}\sum_{k}P_{hk}P_{kj}=\mathbbm{P}\left[\mbox{Transition
to}\,X_{i}\,|\,\mbox{Current Node}=X_{j}\,\text{and}\,3\,\text{steps}\right].$

It follows by induction that

(3) $\displaystyle P_{i,j}^{N}=\mathbbm{P}[\text{Transition to}\,X_{j}\,|$
$\displaystyle\mbox{Current Node}=X_{i}\,\text{and}\,N\,\text{steps}].$

Equation 3 gives us the general formula for the probability of a particle to
travel diffusively from node $i$ to node $j$ in $N$ steps.

To compute the expected value of the average number of steps taken during
diffusion between nodes $i$ and $j$ using Eq. 3, we must specify a new
transition matrix in which there is zero probability of leaving the target
node once the particle has arrived. We denote by
${\color[rgb]{0,0,0}{P_{r}}}\in\mathds{R}^{(n-1)\times(n-1)}$ the new
transition matrix obtained by removing the ${\color[rgb]{0,0,0}{r}}^{th}$ row
and column of the original matrix $P$, where ${\color[rgb]{0,0,0}{r}}$ is our
target node. ${\color[rgb]{0,0,0}{P_{r}}}$ is constructed as follows:

(4)
${\color[rgb]{0,0,0}{P_{r}}}=\begin{pmatrix}\begin{matrix}P_{1:{\color[rgb]{0,0,0}{r}}-1,1:{\color[rgb]{0,0,0}{r}}-1}\end{matrix}&\vline&P_{{\color[rgb]{0,0,0}{r}}+1:n,1:{\color[rgb]{0,0,0}{r}}-1}\\\
\hline\cr
P_{1:{\color[rgb]{0,0,0}{r}}-1,{\color[rgb]{0,0,0}{r}}+1:n}&\vline&\begin{matrix}P_{{\color[rgb]{0,0,0}{r}}+1:n,{\color[rgb]{0,0,0}{r}}+1:n}\end{matrix}\end{pmatrix}.$

From here we may generate an expected value of the amount of time spent on
node $j$ given initial condition $i$ and target node
${\color[rgb]{0,0,0}{r}}$:

(5)
$\displaystyle\left[\;\sum^{\infty}_{k=0}P^{k}_{\color[rgb]{0,0,0}{r}}\right]_{i,j}$
$\displaystyle=[I+{\color[rgb]{0,0,0}{P_{r}}}+{\color[rgb]{0,0,0}{P_{r}}}^{2}+{\color[rgb]{0,0,0}{P_{r}}}^{3}...]_{ij}=[(I-{\color[rgb]{0,0,0}{P_{r}}})^{-1}]_{ij},$
$\displaystyle=\mathbb{E}[\mbox{\\# of visits to $j$ }|\mbox{ Initial
Node}=i]=:\left({\color[rgb]{0,0,0}{N_{r}}}\right)_{i,j},$

where the matrix ${\color[rgb]{0,0,0}{N_{r}}}$ denotes the fundamental matrix
for target condition ${\color[rgb]{0,0,0}{r}}$. The $(i,j)^{\text{th}}$
element of ${\color[rgb]{0,0,0}{N_{r}}}$ represents the average number of
steps a particle will take to node $j$ while diffusing from $i$ to
${\color[rgb]{0,0,0}{r}}$. Note that the formula for
${\color[rgb]{0,0,0}{N_{r}}}$ requires a matrix inversion, for which we will
consider efficient computational strategies later on.

The MFPT is then obtained by summing over row $i$ of
${\color[rgb]{0,0,0}{N_{r}}}$ to find the average amount of time spent
diffusing over the entire network given initial condition $i$ and target
condition ${\color[rgb]{0,0,0}{r}}$:

(6)
${\color[rgb]{0,0,0}{{\color[rgb]{0,0,0}{M_{ir}}}:=({\color[rgb]{0,0,0}{N_{r}}}\cdot\mathbf{Q})_{i}=\sum_{j}\mathbb{E}[\text{\\#
of visits to $j$ }|\text{ Initial Node}=i]\cdot\mathbb{E}[\text{time at
$j$}]}},$

which yields the MFPT matrix $M\in\mathbb{R}^{n\times n}$. The
$(i,r)^{\text{th}}$ entry of $M$ gives the MFPT to go from node $i$ to node
${\color[rgb]{0,0,0}{r}}$.

###  2.3 Mean transit time

The MFPT matrix may be used to assess the average time required to reach a
given node. To compute the mean transit time, or global mean first passage
time (GMFPT) to a node, we average the MFPT to that node with equal weight
over all initial conditions, following [20, 21, 15, 4]. In particular, to
compute the mean transit time to node $r$, we set it as the target condition
as described above. Letting $\widehat{\mbox{M}}_{i}$ denote the $i^{th}$ row
of the MFPT matrix, the GMFPT to node $r$ is given by:

(7) $C_{r}:=\frac{1}{n}\sum_{i}\mathbb{E}[\mbox{travel time to $r$ }|\mbox{
Initial Node}=i]=\frac{1}{n}\cdot\widehat{\mbox{M}}_{r}^{T}\cdot\mathds{1}.$

To provide intuition for this formula, note that the GMFPT is simply the
average of the MFPT of all possible initial conditions to arrive at node $r$.
The $C_{r}$ values at different nodes give a measure of their relative
integration within the network, i.e. smaller values of $C_{r}$ correspond to
nodes that are easily reached from the rest of the network, whereas larger
values of $C_{r}$ correspond to nodes that take longer times to reach.

###  2.4 Edge Importance

In addition to the mean transit time, we may use the MFPT to evaluate the
importance of a given connection between two nodes. To obtain a quantity that
characterizes the overall network, we average $C_{r}$ over all nodes $r$. This
quantity was introduced previously in [8], where it is referred to as the
Target Averaged Global Mean First Passage Time (TAGMFPT). In the notation of
this article, the TAGMFPT is defined as

(8)
$\mbox{TAGMFPT}:=\frac{1}{n}\sum_{r}C_{r}=\frac{1}{n^{2}}\sum_{r}{\color[rgb]{0,0,0}{\widehat{\mbox{M}}_{r}^{T}}}\cdot\mathds{1}.$

We may recompute the TAGMFPT for the same network but after removing the
connection between nodes $i$ and $j$. We shall call this new value
$\mbox{TAGMFPT}^{\star}_{i,j}$, as it is derived using
${\color[rgb]{0,0,0}{N^{\star}_{r}}}$ described in Eq. 13, where
${\color[rgb]{0,0,0}{N^{\star}_{r}}}$ is the given fundamental matrix for
target condition ${\color[rgb]{0,0,0}{r}}$ for our network with one connection
removed. We can look at the difference between these two values to assess the
importance of the connection should it be removed from the network:

(9) $E_{i,j}:=\mbox{TAGMFPT}^{\star}_{i,j}-\mbox{TAGMFPT},$

where $E_{i,j}$ is the importance of the connection between nodes $i$ and $j$.
For intuition, if $E_{i,j}>0$ we can assume that removing the $i,j^{th}$ edge
makes the average amount of steps taken for any path increase.

###  2.5 Synthetic Networks

To benchmark the method, we first applied it to various synthetic networks. We
created simple networks using Delaunay triangulation of a random distribution
of points within a rectangular domain, and then took the Voronoi network (the
dual graph). These networks bear some similarities to biological ER networks
in that the majority of junctions within the network have degree three.
Moreover, the large porous areas created by rings of junctions are reminiscent
of the ER morphology. These properties make them a useful benchmark for our
purposes. In Fig. 1 we show an example of these synthetic networks.

Figure 1: Global mean first massage time (GMFPT) and edge importance of
synthetic networks. (A and B) The GMFPT for the Voronoi network (A) and the
augmented honeycomb network (B). (C and D) Edge importance for the Voronoi
network (C) and the augmented honeycomb network (D). Critical edges in panel
(C) and (D) are labeled in black. The length scale for each network is chosen
so that the mean edge length is 1 µm. Here and throughout the manuscript,
times are given in units of seconds.

We begin by testing our method on networks with a fixed number of connections
per node to quantify how the degree of a given junction affects MFPTs across
the network. The Voronoi network as shown in Fig. 1A is very close to having
all degree-three nodes, except for nodes on the boundary that may have fewer
edges.

As a benchmark comparison similar to that reported in [8, 4], we generated an
augmented, non-symmetric honeycomb structured network through randomly
removing connections from its corresponding symmetric form (see Fig. 1B).
These honeycomb networks have a propensity to form long and isolated branches,
which are significant when testing nodes for GMFPT. The honeycomb networks are
created by first generating the centroid of each hexagon within the network.
The middle row of hexagons is always made odd to preserve symmetry of the
network.

To analyze the list of synthetic networks described above, we make use of the
MPFT described in Section 2.2. For the Voronoi network (Fig. 1A), we
quantified its GMFPT using Eq. 7 and found that the GMFPT of a given node is
most dependent on centrality (distance from the center of network). This is
due to the fact that these networks are mainly composed of degree three nodes
and do not form isolated branches extending from the main network.

The augmented honeycomb network (Fig. 1B) has a feature that distinguishes it
from the other synthetic networks. For this particular type of network, all
features are the same as the original honeycomb except edges were removed with
a probability of 0.2. This means that a node within the interior of the
network has a $0.2$ probability to have two connections, and $0.2^{2}=0.04$
probability to have only one connection. This creates large isolated branches
within the network, which heavily influence the GMFPT per node (see Fig. 1B),
while maintaining overall network properties such as the TAGMFPT.

To further investigate the importance of an edge, we repeat the centrality
analysis. For the Voronoi network we find that the importance of an edge is
mainly influenced by two traits: 1) how central the connection is, and 2) if
removing that connection will increase the minimum path between any two nodes
greatly (see Fig. 1C).

The augmented honeycomb network is unique from the synthetic Voronoi network.
This is because there exist connections within the network such that if they
were to be removed, two disconnected networks would be created. We shall call
these connections critical edges. We can check for these critical edges by
solving for the condition Eq. 16 as described in Section 2.6. In Fig. 1D we
label these edges in black.

###  2.6 The Sherman Morrison formula

To study the effect of temporal changes in the network, such as removal of an
edge, in practice one needs an efficient method to compute the MFPT matrix Eq.
6. The Sherman-Morrison (SM) formula is a well-known formula for the inverse
of a matrix with a rank-1 update. Using this formula we may efficiently invert
perturbations of the MFPT matrix e.g. for swapped target conditions. In
particular, given vectors $u,v\in\mathds{R}^{n}$ and an invertible matrix
$A\in\mathds{R}^{n\times n}$ with $(A+uv^{T})^{-1}$ invertible, the SM formula
yields the inverse of $A+uv^{T}$ as

(10) $(A+uv^{T})^{-1}=A^{-1}-\frac{A^{-1}uv^{T}A^{-1}}{1+v^{T}A^{-1}u}$

Note the assumption that $(A+uv^{T})^{-1}$ is invertible, or equivalently
${1+v^{T}A^{-1}u}\neq 0$ upon taking the determinant. We will discuss the
importance of this criterion later on.

We denote by ${\color[rgb]{0,0,0}{N_{r}}}$ the fundamental matrix of our
original network and ${\color[rgb]{0,0,0}{N^{\star}_{r}}}$ the fundamental
matrix of a network with one edge removed. Using Eq. 5 we can make a series of
rank 1 updates to ${\color[rgb]{0,0,0}{P_{r}}}$ to create
${\color[rgb]{0,0,0}{P^{\star}_{r}}}$, and then using the vectors of the rank
1 updates, we can apply the SM formula to compute the same rank one updates on
the inversion ${\color[rgb]{0,0,0}{N_{r}}}$ to create
${\color[rgb]{0,0,0}{N^{\star}_{r}}}=(I-{\color[rgb]{0,0,0}{P^{\star}_{r}}})^{-1}$.
For a simple change like an edge removal, our transition matrices for each
network, ${\color[rgb]{0,0,0}{P_{r}}}$ and
${\color[rgb]{0,0,0}{P^{\star}_{r}}}$ are very similar. If an edge was removed
between nodes $i$ and $j$, then the difference between
${\color[rgb]{0,0,0}{P_{r}}}$ and ${\color[rgb]{0,0,0}{P^{\star}_{t}}}$ would
only be in the $i^{th}$ and $j^{th}$ rows and columns. We let
$e_{k}\in\mathds{R}^{n-1}$ be the $k^{th}$ standard basis vector and
$p_{i},p^{\star}_{i}\in\mathds{R}^{n-1}$ be the $i^{th}$ column of
${\color[rgb]{0,0,0}{P_{r}}}$ and ${\color[rgb]{0,0,0}{P^{\star}_{r}}}$
respectively. Similarly, we denote by $b_{i},b^{\star}_{i}\in\mathds{R}^{n-1}$
the $i^{th}$ rows of ${\color[rgb]{0,0,0}{P_{r}}}$ and
${\color[rgb]{0,0,0}{P^{\star}_{r}}}$ respectively. Given that the only
difference between ${\color[rgb]{0,0,0}{P_{r}}}$ and
${\color[rgb]{0,0,0}{P^{\star}_{r}}}$ is to the rows and columns associated
with the nodes connected to the edge we wish to remove, we can easily equate
the two. This is accomplished by swapping out the necessary columns and rows:

(11)
${\color[rgb]{0,0,0}{P^{\star}_{r}}}={\color[rgb]{0,0,0}{P_{r}}}+(p^{\star}_{i}-p_{i})e^{T}_{i}+(p^{\star}_{j}-p_{j})e^{T}_{j}+e_{i}(b^{\star}_{i}-b_{i})^{T}+e_{j}(b^{\star}_{j}-b_{j})^{T}$

We subsequently introduce two vectors $\rho_{i}=p^{\star}_{i}-p_{i}$ and
$\beta_{i}=b^{\star}_{i}-b_{i}$ to rewrite (11) in a simpler form:

(12)
${\color[rgb]{0,0,0}{P^{\star}_{t}}}={\color[rgb]{0,0,0}{P_{t}}}+\rho_{i}e^{T}_{i}+\rho_{j}e^{T}_{j}+e_{i}\beta_{i}^{T}+e_{j}\beta_{j}^{T}.$

We note that (12) is now a series of rank 1 updates on
${\color[rgb]{0,0,0}{P_{r}}}$ and the fundamental matrix
${\color[rgb]{0,0,0}{N^{\star}_{r}}}$ may be calculated by

(13)
${\color[rgb]{0,0,0}{N^{\star}_{r}}}=(I-{\color[rgb]{0,0,0}{P_{r}}}-\rho_{i}e^{T}_{i}-\rho_{j}e^{T}_{j}-e_{i}\beta_{i}^{T}-e_{j}\beta_{j}^{T})^{-1}.$

We can then apply the SM formula to extract
${\color[rgb]{0,0,0}{N^{\star}_{r}}}$ from ${\color[rgb]{0,0,0}{N_{r}}}$ in a
computationally efficient manner, using the vectors provided. From here we may
easily compute $M^{\star}$, which is our MFPT matrix for the same network
except with our desired edge removal.

When dealing with inverting a matrix we must also check that the matrix is
indeed invertible. We must check that the determinant of our new matrix
satisfies

(14) $\det(A+uv^{T})=(1+v^{T}A^{-1}u)\det(A).$

Given that ${\color[rgb]{0,0,0}{N^{\star}_{r}}}$ is a series of rank 1 updates
on ${\color[rgb]{0,0,0}{N_{r}}}$ as described in Eq. 13, we may compute the
determinant:

(15)
$\det(I-{\color[rgb]{0,0,0}{P^{\star}_{r}}})=(1+e^{T}_{i}{\color[rgb]{0,0,0}{N_{r}}}^{1}\rho_{i})(1+e^{T}_{j}{\color[rgb]{0,0,0}{N_{r}}}^{2}\rho_{j})(1+\beta_{i}^{T}{\color[rgb]{0,0,0}{N_{r}}}^{3}e_{i})(1+\beta_{j}^{T}{\color[rgb]{0,0,0}{N_{r}}}^{\star}e_{j})\det(I-{\color[rgb]{0,0,0}{P_{r}}}).$

Note that
${\color[rgb]{0,0,0}{N_{r}}}^{1},{\color[rgb]{0,0,0}{N_{r}}}^{2},{\color[rgb]{0,0,0}{N_{r}}}^{3},{\color[rgb]{0,0,0}{N_{r}}}^{\star}$,
define the fundamental matrix after each rank one update to
$I-{\color[rgb]{0,0,0}{P_{r}}}$. We assume that ${\color[rgb]{0,0,0}{N_{r}}}$
is invertible, which is equivalent to $\det(I-{\color[rgb]{0,0,0}{P_{r}}})\neq
0$, meaning that the original network is irreducible. This simplifies the
condition of checking if the new ${\color[rgb]{0,0,0}{N^{\star}_{r}}}$ is
invertible as we only need to show that:

(16)
$(1+e^{T}_{i}{\color[rgb]{0,0,0}{N_{r}}}^{1}\rho_{i})(1+e^{T}_{j}{\color[rgb]{0,0,0}{N_{r}}}^{2}\rho_{j})(1+\beta_{i}^{T}{\color[rgb]{0,0,0}{N_{r}}}^{3}e_{i})(1+\beta_{j}^{T}{\color[rgb]{0,0,0}{N_{r}}}^{\star}e_{j})\neq
0$

Given that we only need to check the rank one updates individually for
invertibility, we may simplify (16) by replacing
${\color[rgb]{0,0,0}{N_{r}}}^{1},{\color[rgb]{0,0,0}{N_{r}}}^{2},{\color[rgb]{0,0,0}{N_{r}}}^{3},{\color[rgb]{0,0,0}{N_{r}}}^{\star}$
with ${\color[rgb]{0,0,0}{N_{r}}}$.

(17)
$(1+e^{T}_{i}{\color[rgb]{0,0,0}{N_{r}}}\rho_{i})(1+e^{T}_{j}{\color[rgb]{0,0,0}{N_{r}}}\rho_{j})(1+\beta_{i}^{T}{\color[rgb]{0,0,0}{N_{r}}}e_{i})(1+\beta_{j}^{T}{\color[rgb]{0,0,0}{N_{r}}}e_{j})\neq
0$

We note that Eq. 17 does not hold in general. If the criterion Eq. 17 is not
satisfied, we deduce that the connection removed to create
${\color[rgb]{0,0,0}{N^{\star}_{r}}}$ results in a reducible network with two
disconnected sub-networks. However, there is one exception: if a connection is
removed that only separates one node from the rest of the network, then the
deduction is false because in the case of a disconnected network Eq. 17 does
not hold. However, we may account for this special case by first testing if
removing a given edge results in a degree zero node. This is computationally
inexpensive, and the combination of checking if Eq. 17 is true and that no
nodes are degree zero is sufficient to conclude that the new network is
irreducible.

The SM formula may also be used to efficiently compute all possible target
conditions of the fundamental matrix. This is explained in detail in Appendix
A.

Figure 2: Computation time and efficiency comparison of the original method
and the optimized method using the SM formula. (A) Computation time of both
methods for the Voronoi network (blue) and the augmented honeycomb network
(red) with varying complexity. The solid lines with triangle markers
correspond to the original method without using the SM formula. The solid
lines with circle markers correspond to the optimized method using the SM
formula. (B) The corresponding speed-up factor for the Voronoi network (blue-
circle) and the augmented honeycomb network (red-triangle) calculated using
data shown in (A).

Using these two applications of the SM formula, we now only need to directly
invert a single matrix, and from this single matrix inverse we can compute any
target condition or simple network change without having to directly invert a
new matrix. To demonstrate the efficiency of this method, we compared the cost
of computing the MFPT matrix using both the original and optimized methods on
networks of increasing size. In Fig. 2A, we show that the computation time
increases significantly using the original method as the number of nodes is
increased for both the Voronoi and the augmented honeycomb networks; whereas
the computation time remains small using the optimized method. To further
quantify the difference between the two methods, we calculated the speed-up
factor for each network complexity. Here the speed-up factor is defined to be
the ratio of the slower to the faster computation time. For both network
types, we see a significant increase in the computation speed using the
optimized method (Fig. 2B) as the network becomes more complex.

##  3 Application to ER networks

To explore the importance of network properties in dynamic biological
organelles, we applied MFPT analysis to a series of ER networks. These
networks expand over an entire cell and consist of hundreds of nodes and
thousands of connections. Instead of attempting to accurately model the
network remodeling dynamics, here we extract the network dynamics directly
from data and compare the observed dynamics to an extremely simple version of
remodeling in which only single edge deletion is allowed. Recall that the
difference in diffusive transport time upon edge deletion is precisely the
quantity measured by the edge importance introduced previously.

Figure 3: Diffusive transport analysis of an S2 cell. (A) The original cell
image and the corresponding (B) extracted network, (C) GMFPT, and (D) Edge
importance. A 5 micron scale bar is plotted at the right bottom corner in
panel (A). The critical edges in panel (D) are shown in black (see also
Appendix C). Figure 4: Centrality analysis of individual node or edge. (A, C,
and E) Centrality of individual node compared to its GMFPT $C_{r}$ for the
Voronoi network with 188 nodes (A), the augmented honeycomb network with 92
nodes (C), and an S2 cell (E). (B, D, and F) Centrality of individual edge
compared to its relative importance for the Voronoi network with 188 nodes
(B), the augmented honeycomb network with 92 nodes (D), and an S2 cell (F).

We used two separate sets of data obtained from experiments on _D.
melanogaster_. First, we examined a series of experimental ER images from an
S2 cell obtained as described in Section B.1. We demonstrate an example of the
ER network extracted from a snapshot (Fig. 3A). Figure 3B shows that the S2
cell has a complicated ER network consisting of different types of regions.
The region around the center of the ER network is dense with nodes highly
connected to each other whereas the region close to the boundary of the
network is relatively sparse. As a result, nodes in the inner region have a
smaller impact on the MFPT compared to those in the outer region (Fig. 3C). We
observe a similar trend in changes in edges in the inner region. That is,
removing edges in the highly connected region would have a minimal effect on
the calculated MFPT. However, in contrast to the nodes, edges with the highest
importance correspond to those that bridge the gap between regions (see the
left bottom corner in Fig. 3D).

Figure 5: Diffusive transport analysis of the NMJ. (A) Original microscopy
image. (B) The corresponding extracted network. (C) GMFPT. (D) Edge
importance. A 2 micron scale bar is plotted at the right bottom corner in
panel (A). In panel (D) the critical edges are plotted in black (see also
Appendix C).

To further quantify the structural properties of the S2 cell ER network, we
performed a centrality analysis, which relates the distance of a node or an
edge from the center of the network to its impact on the MFPT. In comparing
centrality of a node in the S2 cell ER network to the synthetic networks, we
found that the S2 cell ER network (Fig. 4E) is most similar to the Voronoi
network (Fig. 4A), which is designed to maximise the degree of each node. It
is also similar due to the large clusters of porous areas across each network.
For all networks there is a positive relationship between $C_{r}$, GMFPT, and
the distance from the center of the network. In examining the structure of
each network, we see that there are fewer edges connected to nodes on the
periphery of the network (i.e. located far from the center of the network)
than nodes close to the center of the network. As a result, it takes a longer
time to reach nodes closer to the boundary.

Figure 6: The TAGMFPT of the S2 cell as a function of time, with 95$\%$
confidence intervals for each data point calculated by averaging over total
number of nodes.

Repeating the centrality analysis for individual edges for all three types of
networks (Fig. 4B, D, and F), we find that the S2 cell ER network demonstrates
a unique relationship between edge importance and distance from center (Fig.
4F). In contrast to the Voronoi and the augmented honeycomb networks, we see
that the edge importance for the S2 cell ER network is more heterogeneous, in
that edges closer to the center have a lower importance on average. We observe
that the region around the center/nucleus of the ER network is a densely
connected area with short edges; therefore, we hypothesize that removing one
edge would not have a significant effect in the transport process of a
particle. As the distance from the center of the network is increased, we
start to see an increasing number of edges that have high importance. We note
that most of these edges serve as bridges that connect the dense central
region of the network to its sparse boundary region. Removing such an edge
could completely separate the ER network into two disconnected pieces.

The second set of ER network data we examined is taken from a D. melanogaster
neuron at the neuromuscular junction (NMJ), with the experimental protocol
described in Section B.2. We apply the same routine to the NMJ ER network. We
use the NMJ image of Fig. 5A and show the corresponding extracted network in
Fig. 5B. The overall size of the NMJ ER network is significantly smaller and
less dense than the S2 networks. Analyzing the GMFPT we find a similar trend
in relation to the distance to the center as for the S2 cell ER network. We
see from Fig. 5C that nodes that are close to the center have smaller GMFPTs
compared to those close to the boundary.

To investigate the effect of changes over time on diffusive transport in the
ER network, we extracted the ER network from timelapse data from both the S2
cell and NMJ. We repeated the calculations of Section 2.5 frame-by-frame to
compute the changes over time in TAGMFPT for the S2 cell (Fig. 6).

Figure 7: Edge importance over time taken from a timeseries of the NMJ. (A)
1${}^{\text{st}}$ frame, (B) 10${}^{\text{th}}$ frame, (C) 50${}^{\text{th}}$
frame, and (D) 100${}^{\text{th}}$ frame. In all panels the critical edges are
shown in black.

In addition, we analyze the changes in the edge importance over time computed
for the NMJ by comparing the $1^{\text{st}}$, $10^{\text{th}}$,
$50^{\text{th}}$, and $100^{\text{th}}$ frames (Fig. 7). We see that most
edges, except the boundary edges, have a low importance value initially (Fig.
7A and B). In this particular example, as time increases the ER network
becomes more fragmented and edge importance increases as indicated by the
additional black edges that appear in Fig. 7C and D. This is due to the fact
that removing any of those edges would disconnect one or more nodes from the
rest of the network.

##  4 Discussion

In this work we have developed a coarse-grained model via random walks along
with an optimized method to simulate the diffusive transport on networks that
are evolving over time. We describe how the mean first passage time may be
used to quantify the influence of the changing network morphologies on
diffusive transport. Whereas on a static network our approach in calculating
the MFPT is similar to several past studies of transport on spatial networks
[4, 22], we describe an efficient computational method that allows us to study
how the diffusive transport is affected by the evolving network structure.

After benchmarking our method using various synthetic networks, we apply it to
ER network data from S2 cells and the NMJ in _D. melanogaster_. We find that
the GMFPT appears to demonstrate a similar trend in both cases: 1) nodes that
are closer to the center tend to be more connected and therefore have a lower
GMFPT 2) nodes that are closer to the boundary tend to have a higher GMFPT as
they are much harder to reach by diffusion. The edge importance reveals more
information about the effect of the overall network properties of the ER on
diffusion. The S2 cell network has two regions with different node densities
that are bridged by relatively few edges. As a result, these edges have a
higher edge importance as deleting one will separate the whole network into
disconnected components with some inaccessible to the particle.

In reviewing the time evolution of the centrality of the NMJ networks, we
report an interesting relation between edges with negative edge importance and
their locations. Figure 7 shows that these edges tend to form small loops and
are connected to critical edges with high edge importance. The observation of
negative edge importance is an instance of Braess’ paradox [3], the
counterintuitive finding that decreasing network capacity can under certain
conditions lead to faster average transit times. Braess’ paradox has been
observed to occur in road networks [18] as well as physics systems such as
electronic circuits [5], and given that similar modeling assumptions (e.g.
Kirchoff’s laws) are used in the present study, it is not altogether
surprising that it may occur in the present context as well. However, studying
the prevalence with which Braess’ paradox occurs in biological networks and
analyzing its consequences for transport efficiency would be an interesting
avenue for further study.

Whereas the GMFPT’s discussed above are obtained by averaging with equal
weight over the MFPT, one could alternatively average with respect to the
stationary distribution (Fig. 8). This yields the expected MFPT for initial
conditions chosen according to the stationary distribution, i.e. the steady-
state distribution for particles diffusing over the network (Fig. 9). As shown
in Fig. 10, the qualitative features of the edge importance are not
significantly changed by this choice of weighting.

We remark that although here we apply the edge importance to biological
networks that change over time, the method is not inherently tied to dynamics
and may be applied to other contexts as well, e.g. to identify critical edges
within transport networks. This is because the principle of using the Sherman-
Morrison formula to efficiently compute the edge importance by removing one
edge at a time applies more generally and is not limited to time-dependent
networks. In addition, because the edge importance is computed independently
for each edge, the algorithm may be parallelized in a relatively
straightforward manner to further accelerate the computation.

There are a number of different directions that this research could take,
including introducing drift to the transport process. Recent experiments have
revealed that transport of particles on ER networks is affected by microtubule
dynamics and is not purely diffusive [14]. Generalizing the methods presented
here to non-diffusive transport would be a promising avenue for future
research.

Figure 8: Stationary distribution of the MFPT for each node in ER networks.
(A) 1${}^{\text{st}}$ frame of the S2 cell timeseries. (B) 1${}^{\text{st}}$
frame of the NMJ timeseries. Figure 9: GMFPT of ER networks obtained using a
nonuniform stationary distribution weighting. (A) 1${}^{\text{st}}$ frame of
the S2 cell timeseries. (B) 1${}^{\text{st}}$ frame of the NMJ timeseries.
Each node is weighted by its corresponding stationary distribution shown in
Fig. 8. Figure 10: Comparison of edge importance of ER networks using
different weightings. (A) 1${}^{\text{st}}$ frame of the S2 cell timeseries
using the uniform weighting. (B) 1${}^{\text{st}}$ frame the S2 cell
timeseries using the nonuniform stationary weighting. (C) 1${}^{\text{st}}$
frame of the NMJ timeseries using the uniform weighting. (D) 1${}^{\text{st}}$
frame of the NMJ timeseries using the nonuniform stationary weighting. The
critical edges are plotted in black.

##  Appendix A Sherman-Morrison formula

The SM formula may be used to compute the fundamental matrices for different
target conditions. That is to say we can go from ${\color[rgb]{0,0,0}{N_{r}}}$
to $N_{i}\;\;\forall i\in\\{1,2,3,...,n\\}$. If we choose
${\color[rgb]{0,0,0}{r}}=n$ we can use ${\color[rgb]{0,0,0}{N_{r}}}$ as our
base matrix upon which we can switch to any new target through a series of the
SM formula applications. This process is slightly more complicated than the
method for computing ${\color[rgb]{0,0,0}{N^{\star}_{r}}}$ because
${\color[rgb]{0,0,0}{P_{r}}}\in\mathds{R}^{(n-1)\times(n-1)}$ requires an
extra row and column in order to swap in the probabilities of the new target
condition. We define
${\color[rgb]{0,0,0}{\widetilde{P_{r}}}}\in\mathds{R}^{n\times n}$ to be
identical to the original $P$ except the ${\color[rgb]{0,0,0}{r}}^{th}$ row
and column will be all zeros, but
$({\color[rgb]{0,0,0}{\widetilde{P_{r}})_{r,r}}}=1$. For
${\color[rgb]{0,0,0}{r}}=n$ we have:

(18)
$\widetilde{P_{n}}=\begin{pmatrix}\begin{matrix}P_{n}\end{matrix}&\vline&\mathbf{0}\\\
\hline\cr\mathbf{0}&\vline&\begin{matrix}1\end{matrix}\end{pmatrix}.$

For $r\neq n$, the rows and columns removed will split the
${\color[rgb]{0,0,0}{P_{r}}}$ matrix into four block matrices along the
$r^{th}$ row and column as follows:

(19)
$\widetilde{P_{r}}=\begin{pmatrix}P^{UL}_{r}&\vline&\mathbf{0}&\vline&P^{UR}_{r}\\\
\hline\cr\mathbf{0}&\vline&1&\vline&\mathbf{0}\\\ \hline\cr
P^{LL}_{r}&\vline&\mathbf{0}&\vline&P^{LR}_{r}\\\ \end{pmatrix},$

where $P^{UL}_{r},P^{LR}_{r}\in\mathbb{R}^{(r-1)\times(r-1)}$ and
$P^{UR}_{r},P^{LL}_{r}\in\mathbb{R}^{(n-r)\times(n-r)}$. We let
$P_{r}\in\mathbb{R}^{(n-1)\times(n-1)}$ denote the original probability matrix
without the r${}^{\text{th}}$ row and column and is defined as

(20)
$P_{r}=\begin{pmatrix}\begin{matrix}P^{UL}_{r}\end{matrix}&\vline&P^{UR}_{r}\\\
\hline\cr
P^{LL}_{r}&\vline&\begin{matrix}P^{LR}_{r}\end{matrix}\end{pmatrix}.$

Having defined the ${\color[rgb]{0,0,0}{\widetilde{P_{r}}}}$ matrices, we may
form rank-1 updates of $\widetilde{P_{n}}$ to obtain
$\widetilde{P_{r}}\;\;\forall r\in\\{1,2,3...,n\\}$. We define vectors $u_{i}$
and $v_{i}$ to be the $i^{\text{th}}$ row and column of our original $P$
except for the $i^{\text{th}}$ element:
$(u_{i})_{i},(v_{i})_{i}=\frac{1}{2}(1-P_{i,i})$, this is to account for the 1
in the zero rows. The updates go as follows:

(21)
$\displaystyle\widetilde{P_{r}}=\widetilde{P_{n}}-e_{r}u^{T}_{r}-v_{r}e^{T}_{r}+e_{n}u^{T}_{n}+v_{n}e^{T}_{n}$

From here we may define $\widetilde{N_{r}}$, which is the original fundamental
matrix except split across the same row and vector as $\widetilde{P_{r}}$:

(22)
$\widetilde{N_{r}}=\begin{pmatrix}N^{UL}_{r}&\vline&\mathbf{0}&\vline&N^{UR}_{r}\\\
\hline\cr\mathbf{0}&\vline&1&\vline&\mathbf{0}\\\ \hline\cr
N^{LL}_{r}&\vline&\mathbf{0}&\vline&N^{LR}_{r}\\\
\end{pmatrix}=\begin{pmatrix}(I-P_{r})^{-1}_{UL}&\vline&\mathbf{0}&\vline&(I-P_{r})^{-1}_{UR}\\\
\hline\cr\mathbf{0}&\vline&1&\vline&\mathbf{0}\\\
\hline\cr(I-P_{r})^{-1}_{LL}&\vline&\mathbf{0}&\vline&(I-P_{r})^{-1}_{LR}\\\
\end{pmatrix},$

and from this matrix one may derive $N_{r}$ through the same process used to
obtain $P_{r}$ from $\widetilde{P_{r}}$.

We may now compute $\widetilde{N_{r}}$ as a series of rank 1 updates on
$\widetilde{P_{n}}$ given by the vectors provided above:

(23)
$\widetilde{N_{r}}=(I-(\widetilde{P_{n}}-e_{r}u^{T}_{r}-v_{r}e^{T}_{r}+e_{n}u^{T}_{n}+v_{n}e^{T}_{n}))^{-1}.$

Using $\widetilde{N_{n}}=(I-\widetilde{P_{n}})^{-1}$ as a base inversion, we
may now apply the 4 rank 1 updates through the SM formula to compute
$\widetilde{N_{r}}$, and by removing the $r^{th}$ row and column we get the
matrix $N_{r}$ from which we can derive the $r^{th}$ row of our matrix $M$.
Iterating over all $r\leq n$ we can derive the entire MFPT matrix.

It is important to note that the order of the rank 1 updates matters. If one
first adds $e_{n}u^{T}_{n}+v_{n}e^{T}_{n}$ before subtracting
$e_{r}u^{T}_{r}+v_{r}e^{T}_{r}$, it will create a non-invertible matrix
halfway through the process and not be able to generate $\widetilde{N}_{r}$.
The intuition behind this is that you are removing the old target condition
before adding the new one, and thus with no target condition you cannot have
mean first passage times, resulting in a non-invertible matrix.

##  Appendix B Experimental materials and methods

###  B.1 S2 Cell Imaging

To generate pQUAST-BiP sfGFP HDEL, the sequence for BiP sfGFP HDEL was
obtained by Dr. James McNew at Rice University and was cloned into pQUAST-attB
(Addgene$\\#$ 104880) by Vector Builder (Chicago, IL).

S2 cells were maintained at $28$ degree Celsius in Schneider’s Drosophila
Medium (Thermo Fisher) with 10$\%$ heat inactivated fetal bovine serum and
$50$ units of penicillin-streptomycin per mL. Cells were transfected using
Effectene (Qiagen, Hilden, Germany) following the manufacturer’s instructions,
with pAC-7-QFBDAD (Addgene$\\#$ 46096) to drive expression of QUAS-BiP sfGFP
HDEL for 24 hours. Cells were transferred to cover glass chambers (Thermo
Scientific, Waltham MA), coated with $0.5$ mg/ml Concanavalin A and left to
settle for 1 hour at room temperature. Cells were then imaged at room
temperature with a 63X (NA1.4) oil immersion objective, using an Airyscan LSM
880 microscope in super-resolution mode and Zen Black software. To record ER
dynamics of _Drosophila_ S2 cells, 100 frames of a single slice were obtained
at a frame rate of 0.65 seconds. For analysis, images were corrected for
photobleaching using the histogram matching method with the imaging software
FIJI [13].

###  B.2 _Drosophila_ larval neuromuscular junction

_Drosophila melanogaster_ were cultured on standard medium at 25 degree
Celsius. To label neuronal ER, Vglut-Gal4 (Bloomington Drosophila Stock Center
(BDSC) $\\#$24635) flies were crossed with UAS BiP sfGFP HDEL(BDSC $\\#$64748)
flies. 3rd instar wandering larvae were dissected in Ca 2+ -free HL3. 1 [6]
and axons were severed from the central nervous system. Larvae were dissected
in glass slides with Press To Seal Silicone Isolators (Grace Bio-labs;
CQS-13R-2.0) filled at their centers with Krayden Dow Sylgard 184 Silicone
(Thermo Fisher Scientific). Small metal pins cut with nail clippers were used
to stretch the larvae and the pins were pressed into the cured Sylgard 184
Silicone. Finally, a coverslip was placed on top of the dissected larvae and
excess HL3. 1 was removed before imaging. Larvae were then imaged at room
temperature with a 63X (NA1.4) oil immersion objective, using an Airyscan LSM
880 microscope in super-resolution mode and Zen Black software. To record ER
dynamics of _Drosophila_ nerve terminals, imaging was performed by taking
Z-stacks of 8 slices with 0.5 micron spacing for 100 frames at a frame rate of
0.73 seconds. Next, images were corrected for photobleaching and muscle
contraction with FIJI using histogram matching and the “StackReg” plugin [16],
respectively.

Finally, maximum projections of images were generated for analysis. The nodal
and edge information are extracted using MATLAB. To extract the network
information for individual frames, we begin by denoising the individual image
through gridded interpolation followed by a background subtraction and guided
filtering. An anisotropic diffusion filter and a threshold adjustment are
applied to the denoised image to enhance the tubular structure. To obtain the
corresponding image skeleton, we make use of the MATLAB functions bwmorph and
bwconncomp to extract the maximum connected component of the skeleton. The
resulting binary skeleton is converted to a graph to obtain the nodal and edge
information for each image.

##  Appendix C Critical edges of ER networks

Critical edges, defined as the edges whose removal would disconnect the
network into multiple components, are shown for two representative ER networks
in Fig. 11.

Figure 11: Critical edges of ER networks. (A) 1${}^{\text{st}}$ frame of the
S2 cell timeseries. (B) 1${}^{\text{st}}$ frame of the NMJ timeseries. The
critical edges are shown in red.

## Acknowledgments

We acknowledge helpful discussions with Steven Petteruti.

## References

  * [1] D. Ando, N. Korabel, K. C. Huang, and A. Gopinathan, Cytoskeletal network morphology regulates intracellular transport dynamics, Biophysical Journal, 109 (2015), pp. 1574–1582.
  * [2] O. Bénichou, C. Chevalier, B. Meyer, and R. Voituriez, Facilitated diffusion of proteins on chromatin, Phys. Rev. Lett., 106 (2011), p. 038102.
  * [3] D. Braess, Über ein paradoxon aus der verkehrsplanung, Unternehmensforschung, 12 (1968), pp. 258–268.
  * [4] A. I. Brown, L. M. Westrate, and E. F. Koslover, Impact of global structure on diffusive exploration of organelle networks, Sci. Rep., 10 (2020).
  * [5] J. E. Cohen and P. Horowitz, Paradoxical behaviour of mechanical and electrical networks, Nature, 352 (1991), pp. 699–701.
  * [6] Y. Feng, A. Ueda, and C.-F. Wu, A modified minimal hemolymph-like solution, hl3. 1, for physiological recordings at the neuromuscular junctions of normal and mutant drosophila larvae, Journal of Neurogenetics, 18 (2004), pp. 377–402.
  * [7] E. Koslover, M. Díaz, D. La Rosa, and A. Spakowitz, Theoretical and computational modeling of target-site search kinetics in vitro and in vivo, Biophysical Journal, 101 (2011), pp. 856–865.
  * [8] E. F. Koslover and A. J. Spakowitz, Force fluctuations impact kinetics of biomolecular systems, Phys. Rev. E, 86 (2012).
  * [9] N. Masuda, M. Porter, and R. Lambiotte, Random walks and diffusion on network, Phys. Rep., 716-717 (2017), pp. 1–58.
  * [10] R. Milo and R. Phillips, Cell biology by the numbers, Garland Science, 2015.
  * [11] L. Niu, T. Ma, F. Yang, B. Yan, X. Tang, H. Yin, Q. Wu, Y. Huang, Z.-P. Yao, J. Wang, Y. Guo, and J. Hu, Atlastin-mediated membrane tethering is critical for cargo mobility and exit from the endoplasmic reticulum, PNAS, 116 (2019), pp. 14029–14038.
  * [12] S. Petteruti, Fast simulation of random walks on the endoplasmic reticulum by reduction to the gambler’s ruin problem, Undergraduate Honors Thesis, (2017).
  * [13] J. Schindelin, I. Arganda-Carreras, E. Frise, V. Kaynig, M. Longair, T. Pietzsch, S. Preibisch, C. Rueden, S. Saalfeld, B. Schmid, et al., Fiji: an open-source platform for biological-image analysis, Nature Methods, 9 (2012), pp. 676–682.
  * [14] L. Stadler, K. Speckner, and M. Weiss, Diffusion of exit sites on the endoplasmic reticulum: A random walk on a shivering backbone, Biophysical Journal, 115 (2018), p. 1552.
  * [15] V. Tejedor, O. Bénichou, and R. Voituriez, Global mean first-passage times of random walks on complex networks, Physical Review E, 80 (2009), p. 065104.
  * [16] P. Thevenaz, U. E. Ruttimann, and M. Unser, A pyramid approach to subpixel registration based on intensity, IEEE transactions on image processing, 7 (1998), pp. 27–41.
  * [17] M. Viana, A. Brown, I. Mueller, C. Goul, E. Koslover, and S. Rafelski, Mitochondrial fission and fusion dynamics generate efficient, robust, and evenly distributed network topologies in budding yeast cells, Cell Syst., 10 (2020), pp. 287–297.
  * [18] H. Youn, M. T. Gastner, and H. Jeong, Price of anarchy in transportation networks: efficiency and optimality control, Physical review letters, 101 (2008), p. 128701.
  * [19] B. Zani and E. Edelman, Cellular bridges: routes for intercellular communication and cell migration, Commun. Integr. Biol., 3 (2010), pp. 215–220.
  * [20] Z. Zhang, J. Guan, W. Xie, Y. Qi, and S. Zhou, Random walks on the apollonian network with a single trap, EPL (Europhysics Letters), 86 (2009), p. 10006.
  * [21] Z. Zhang, Y. Zhang, S. Zhou, M. Yin, and J. Guan, Influences of degree inhomogeneity on average path length and random walks in disassortative scale-free networks, Journal of mathematical physics, 50 (2009), p. 033514.
  * [22] C. Zubenelgenubi, A. Brown, S. Mogre, L. Westrate, and E. Koslover, Diffusive search and trajectories on tubular networks: a propagator approach, The European Physical Journal E, 44 (2021), p. 80.