Proportional mean model for panel count data with multiple modes of recurrence
Sreedevi E. P.111 Corresponding Author Email<EMAIL_ADDRESS>and
Sankaran P. G.2

1 SNGS College, Pattambi.

2 Cochin University of Science and Technology, Cochin.

Abstract: Panel count data is common when the study subjects are exposed to
recurrent events, observed only at discrete time points. In this article, we
consider the regression analysis of panel count data with multiple modes of
recurrence. We propose a proportional mean model to estimate the effect of
covariates on the underlying counting process due to different modes of
recurrence. The simultaneous estimation of baseline cumulative mean functions
and regression parameters of $(k>1)$ recurrence modes are studied in detail.
Asymptotic properties of the proposed estimators are also established. A Monte
Carlo simulation study is carried out to validate the finite sample behaviour
of the proposed estimators. The methods are applied to a real data arising
from skin cancer chemoprevention trial.

Key words : Counting process, panel count data, proportional mean model,
pseudo likelihood, recurrent events.

##  1\. Introduction

In many longitudinal studies on recurrent events in lifetime data analysis,
instead of observing the time to occurrence of event, we may only observe the
number of events experienced by a subject in a given period of time. If each
subject can be observed at more than one time points, the number of events
between two successive observation times is available. The data obtained in
this form is known as panel count data (Kalbfleisch and Lawless,1985; Sun,
2009). Panel count data frequently arise in many fields such as clinical
trials, epidemiological studies and engineering, when continuous follow-up to
obtain exact event times of each subject is infeasible or too costly (Chiou et
al., 2019). Some authors refer panel count data as interval count data or
interval censored recurrent event data (Lawless and Zhan, 1998; Thall, 1988).
An extensive review of panel count data is given in Sun and Zhao (2013). Note
that when the subjects can be observed only at a single time point, we obtain
current status data which is explored in the monograph by Sun (2007).

The standard methods in the analysis of panel count data are focused on the
mean function or the rate function of the underlying recurrent event process.
An estimator for the mean function based on isotonic regression theory is
developed by Sun and Kalbfleisch (1995). Wellner and Zhang (2000) discussed
likelihood based nonparametric estimation methods for the mean function and
proposed a nonparametric maximum likelihood estimator (NPMLE) and a
nonparametric maximum pseudo likelihood estimator (NPMPLE) for the same.
Wellner and Zhang (2000) also studied the asymptotic properties of both NPMPLE
and NPMLE. Thall and Lachin (1988) and Lawless and Zhan (1998) considered the
analysis of panel count data using rate functions. Some of the recent
developments in the analysis of panel count data include Xu et al. (2017) and
Chiou et al. (2018) among others.

In panel count data, it is common to observe a covaraite vector $Z$ for each
subject which affect the underlying counting process of recurrent events. Two
different approaches employed for the analysis of regression models for panel
count data are either by using maximum likelihood methods or by applying the
generalized estimating equation approach. Some important developments in this
area include Sun and Wei (2000), Wellner and Zhang (2007), Zhang (2002) and Hu
et al. (2003). Regression analysis of panel count data with informative
observation times is considered by Haung et al. (2006), Sun et al. (2007) and
Zhao and Tong (2011). Covaraites with measurement error for panel count data
was studied by Kim (2007). Recently Chiou et al. (2019) reviewed various
semiparametric regression modelling approaches for panel count data using R
programming language.

When the study subjects are exposed to recurrent events of several types, we
observe the recurrence due to each possible mode (cause) of recurrence at
different observation times. As a result, we obtain panel count data with
multiple failure modes. For example consider the data on skin cancer
chemoprevention trial discussed in Sun and Zhao (2013). The cancer recurrences
of 290 patients with a history of non-melanoma skin cancers are observed at
different monitoring times. The types of cancers are classified into Basal
cell carcinoma and Squamous cell carcinoma and the recurrences due to both
types of cancers at each monitoring time are observed for each individual.
Covariate information on age, gender, number of prior tumours and DFMO status
is also observed for each individual. Accordingly, we have panel count data
with multiple modes of recurrence with covariates. A detailed analysis of the
data is given in Section 5.

Even though panel count data was a topic of research interest from last two
decades, only a sparse amount of literature is available on panel count data
with multiple failure modes. Sreedevi and Sankaran (2020) developed an
estimator for cause specific mean function and Sankaran et al. (2021) studied
cause specific rate functions of the underlying recurrent event processes when
subjects are exposed to more than one recurrence mode. Both these works
considered data without covariates. Regression analysis of panel count data
with multiple failure modes is not studied yet. Motivated by this, in this
article we propose a proportional mean model to estimate the regression
parameters and baseline cumulative mean functions of panel count data exposed
to more than one mode of recurrence.

The rest of the article is organized as follows. In Section 2, we propose a
new proportional mean model to estimate the baseline cumulative mean functions
and regression parameters due to each mode of recurrence simultaneously . A
simple iterative algorithm is derived for the estimation. Asymptotic
properties of the proposed estimators are established in Section 3. In Section
4, the finite sample behaviour of the proposed estimators is validated through
a Monte Carlo simulation study. The proposed procedures are illustrated using
a real data on skin cancer chemoprevention trial in Section 5. Finally,
Section 6 gives concluding remarks with a discussion on possible future works.

##  2\.  The proportional mean model

Consider a study on $n$ individuals exposed to the recurrent events due to
$\\{1,2,...,k\\}$ different causes. Assume that the event process is observed
only at a sequence of random monitoring times. Consequently, the counts of the
event recurrences due to each mode in between the observation times are only
available; the exact recurrence times remain unknown.

Define a counting process $N_{j}(t)=\\{N_{j}(t);t\geq 0\\}$ where $N_{j}(t)$
denote the number of recurrences of the event due to cause $j$ upto time $t$.
Now, $E(N_{j}(t))=\Lambda_{j}(t)$ for $j=1,2,..,k$ denote the expected number
of cumulative events due to cause (mode) $j$ upto time $t$. The function
$\Lambda_{j}(t)$ is the mean function of the counting process $N_{j}(t)$ and
can be termed as the cause specific mean functions (Sreedevi and Sankaran,
2020). Assume that, corresponding to each subject we observe a $d\times 1$
vector of covariates denoted by $Z$. Our interest is to study
$E(N_{j}(t)|Z)=\Lambda_{j}(t|Z)$ for $j=1,2,..,k$, the expected number of
cumulative events due to cause $j$ upto time $t$ conditionally on covariate
vector $Z$. To estimate the effect of covariate vector $Z$ on lifetime $T$, we
propose the proportional mean model given by

$\Lambda_{j}(t|Z)=\Lambda_{0j}(t)\text{exp}(\beta_{j}^{{}^{\prime}}Z)~{}~{}~{}~{}~{}j=1,2,...,k$
(1)

where $\Lambda_{0j}(.)$ is the completely unspecified baseline mean function
and $\beta_{j}$ is the $d\times 1$ vector of regression parameters
corresponding to cause $j$. When $k=1$, the model in Eq. (1) reduces to the
proportional mean model for panel count data with a single mode of recurrence
studied by Sun and Wei (2000), Zhang (2002) and Wellner and Zhang (2007).

Now, we discuss the structure of panel count data with covariates which is
exposed to multiple modes of recurrence. Let $M$ be an integer valued random
variable denoting the number of observation times which may be different for
each individual and $\underline{T}=\\{T_{m,p},p=1,2,...,m;m=1,2,...\\}$ be the
set of observation times. Now $T_{m,p-1}\leq T_{m,p}$ for $p=1,2,...,m$ and
for all possible values of $m$. Assume that $N_{j}(t)$ and $(M,\underline{T})$
are independent. Let $N_{M,p}^{j}$ denote the number of recurrences of the
event due to cause $j$ upto monitoring time $M$ for $p=1,2,..,M$ and
$j=1,2,...,k$. with $N_{M,p}^{j}=N_{j}(T_{M,p})$. For each subject we also
observe a $d\times 1$ vector of covariates $Z$. Now we observe $n$ i.i.d.
(independent and identically distributed) copies of
$\\{M,T_{M,p},N_{M,p}^{1},...,N_{M,p}^{k},Z\\}$, $p=1,2,...,M$. Accordingly,
observed data will be of the form
$\\{m_{i},t_{m_{i},p},n_{m_{i},p}^{1},...,n_{m_{i},p}^{k},z_{i}\\}$,
$p=1,2,...,m_{i}$ and $i=1,2,...,n$.

The regression analysis of panel count data based on maximum likelihood
methods with a single failure mode is explored by Zhang (2002) and Wellner and
Zhang (2007). Sreedevi and Sankaran (2020) studied panel count data with
multiple modes of recurrence and developed a pseudo likelihood function for
the observed data and derived an isotonic regression estimator (IRE) for cause
specific mean functions. They constructed a pseudo likelihood function for the
observed data (which does not involve covariates) by assuming that the
successive counts of the recurrent event process $N_{j}(t)$ are independent
random variables. We extend their derivation of pseudo likelihood function
into a scenario with covariates. We estimate both $\Lambda_{0j}(.)$ and
$\beta_{j}$ simultaneously as the values that maximize the pseudo likelihood.
Under the assumption that the underlying counting process $N_{j}(t)$ is a non-
homogeneous Poisson process with conditional mean function given in Eq. (1),
for $j=1,2,..,k$, we obtain

$P(N_{j}(t)=m|Z)=\frac{(\Lambda_{0j}\text{exp}(\beta_{j}^{{}^{\prime}}Z))^{m}\text{exp}(\Lambda_{0j}\text{exp}(\beta_{j}^{{}^{\prime}}Z))}{m!}~{}~{}~{}\text{for}~{}~{}~{}m=0,1,2,....$
(2)

When the $k>1$ modes of recurrence are independent, the pseudo log likelihood
function of the observed data can be written as

$l_{n}(\beta_{j},\underline{\Lambda}_{0j},\underline{X})=\sum\limits_{j=1}^{k}l_{nj}(\beta_{j},\Lambda_{0j},\underline{X}),$
(3)

where $\underline{X}$ is the observed data given by
$\underline{X}=\\{M,T_{M,p},N_{M,p}^{1},...N_{M,p}^{k},Z\\}$ and
$l_{nj}(\Lambda_{0j},\underline{X})$ is the log likelihood corresponding to
$j$ th cause. By extending the results in Zhang (2002) for panel count data
with single mode of failure, after ignoring the insignificant parts in the
estimation of $\beta_{j}$’s and $\Lambda_{0j}$’s,
$\l_{nj}(\beta_{j},\Lambda_{0j},\underline{X})$ is given by

$l_{nj}(\beta_{j},\Lambda_{0j},\underline{X})=\sum\limits_{i=1}^{n}\sum\limits_{p=1}^{M_{i}}[N^{j}_{M_{i},p}\text{log}\Lambda_{0j}(T_{M_{i},p})+N^{j}_{M_{i},p}(\beta_{j}^{{}^{\prime}}Z_{i})-\Lambda_{0j}(T_{M_{i},p})\text{exp}(\beta_{j}^{{}^{\prime}}Z_{i})]~{}~{}~{}\text{for}~{}~{}~{}j=1,2,...,k$
(4)

where $M_{i}$ is the number of observation times,
$T_{M_{i},p},~{}p=1,2,..,M_{i}$, the different observation times and
$N^{j}_{M_{i},p},p=1,2,...,M_{i},j=1,2,..,k$ the number of recurrences of the
event due to cause $j$ for $i$ th individual. We assume that given the
covariate vector $Z$, the the distributions of $T$ and $M$ are independent of
$\beta_{j}$ and $\Lambda_{0j}$. We maximize the log pseudo likelihood given in
Eq. (4) to obtain the estimators of $\beta_{j}$ and $\Lambda_{0j}$.
Now we can discuss the computational procedures. Based on the observed data
$\underline{X}$ discussed above, we define the following terms. Let $I(A)$ be
the indicator function of the set $A$ and $s_{1}<s_{2}<...<s_{r}$ be the
distinct ordered observation time points in the set
$\\{T_{M_{i},p},p=1,2,...,M_{i},i=1,2,...,n$}. For $q\in\\{1,2,...,r\\}$ and
for any particular cause of recurrence $J$, define

$b_{qj}=\sum\limits_{i=1}^{n}\sum\limits_{p=1}^{M_{i}}I[{T_{M_{i},p}}=s_{q};J=j]$
(5)

the number of observations made at $s_{q}$ due to cause $j$ and

$\bar{n}_{qj}=\frac{1}{b_{qj}}\
\sum\limits_{i=1}^{n}\sum\limits_{p=1}^{M_{i}}N^{j}_{M_{i},p}I[{T_{M_{i},p}}=s_{q};J=j]$
(6)

as the mean value of the recurrences made at $s_{q}$ due to cause $j$ for
$j=1,2,...,k$. Also define

$V_{qj}(\beta_{j},Z)=\frac{1}{b_{qj}}\sum\limits_{i=1}^{n}\sum\limits_{p=1}^{M_{i}}\text{exp}(\beta_{j}^{{}^{\prime}}Z_{i})I[{T_{M_{i},p}}=s_{q};J=j]$
(7)

and

$W_{qj}(\beta_{j},Z,N)=\frac{1}{b_{qj}}\sum\limits_{i=1}^{n}\sum\limits_{p=1}^{M_{i}}{T_{M_{i},p}}(\beta_{j}^{{}^{\prime}}Z_{i})I[{T_{M_{i},p}}=s_{q};J=j].$
(8)

Now we can rewrite the log pseudo likelihood for $j$ th mode of recurrence
given in Eq.(4) as

$l_{nj}(\beta_{j},\Lambda_{0j}|\underline{X})=\sum\limits_{q=1}^{r}b_{qj}[\bar{n}_{qj}\text{log}\Lambda_{0j}(s_{q})-V_{qj}(\beta_{j},Z)\Lambda_{0j}(s_{q})+W_{qj}(\beta_{j},Z,N_{j})].$
(9)

We maximize Eq.(9) to obtain the estimates of $\beta_{j}$ and $\Lambda_{j}(.)$
for $j=1,2,...,k$. The obtained semiparametric maximum pseudo likelihood
estimators will be the values of parameters that maximize (9) over the set
$R^{d}\times\Omega^{+}$ where $R$ is the set of real numbers and
$\Omega^{+}=\\{(y_{1},y_{2},...,y_{r})\in R^{d}:y_{1}\leq y_{2}\leq...\leq
y_{r}\\}$. The estimators can be obtained as

$(\widehat{\beta_{j}},\widehat{\Lambda_{0j}})=\underset{(\beta_{j},\Lambda_{0j})\in
R^{d}\times\Omega^{+}}{\mathrm{argmax}}{l_{nj}(\beta_{j},\Lambda_{0j}|\underline{X})}.$
(10)

To solve the optimisation problem numerically, we first choose an initial
value of $\beta_{j}$, say $\beta_{j}^{0}$. Now for a fixed $\beta_{j}$, the
estimator of $\Lambda_{0j}$ can be obtained as
$\hat{\Lambda}_{0j}(\beta_{j}^{0})=\underset{\Lambda_{0j}\in\Omega^{+}}{\mathrm{argmax}}(l_{nj}^{*}(\Lambda_{0}j|\beta^{0}_{j},\underline{X}))$
where

$l_{nj}^{*}(\Lambda_{0}j|\beta_{j},\underline{X})=\sum\limits_{q=1}^{r}b_{qj}[\bar{n}_{qj}\text{log}\Lambda_{0j}(s_{q})-V_{qj}(\beta_{j},Z)\Lambda_{0j}(s_{q})].$
(11)

Let $\Lambda_{0j}^{0}$ be the solution of Eq. (11). Now using the estimated
value of $\Lambda_{0j}^{0}$, we can find the updated estimate of $\beta_{j}$
as $\hat{\beta}_{j}(\Lambda_{0j}^{0})=\underset{\beta_{j}\in
R^{d}}{\mathrm{argmax}}(l_{nj}^{**}(\beta_{j}|\Lambda_{0j}^{0},\underline{X})$
where

$l_{nj}^{**}(\beta_{j}|\Lambda_{0}j,\underline{X})=\sum\limits_{q=1}^{r}b_{qj}[W_{qj}(\beta_{j},Z,N_{j})-V_{qj}(\beta_{j},Z)\Lambda_{0j}(s_{q})].$
(12)

The process is continued until the estimators converge. The convergence
criteria can be chosen as

$\displaystyle\left\lvert\frac{l_{nj}^{(h+1)}-l_{nj}^{(h)}}{l_{nj}^{(h)}}\right\rvert\leq\epsilon,$
(13)

where $l_{nj}^{(h)}=l_{n}(\beta_{j}^{(h)},\Lambda_{0j}^{(h)})$ for
$h=0,1,2,....$.
To estimate $(\beta_{j},\Lambda_{0j})$ for $j=1,2,...,k$, the computational
algorithm can be summarised as follows

  1. Step 1:

Choose an initial value $\beta_{j}$ say $\beta_{j}^{0}$.

  2. Step 2:

For the given $\beta_{j}^{h}$, compute $\Lambda_{0j}^{h}$ as the maximum
argument of Eq. (11), given by

$\hat{\Lambda}_{0j}^{h}(\beta_{j}^{h})=\underset{\Lambda_{0j}\in\Omega^{+}}{\mathrm{argmax}}(l_{nj}^{*}(\Lambda_{0}j|\beta^{h}_{j},\underline{X}).$

  3. Step 3:

Update the estimate of $\beta_{j}^{h}$, using the estimate of
$\Lambda_{0j}^{h}$ obtained in Step 2, as the maximum argument of Eq. (12),
given by

$\hat{\beta}_{j}(\Lambda_{0j}^{h})=\underset{\beta_{j}\in
R^{d}}{\mathrm{argmax}}(l_{nj}^{**}(\beta_{j}|\Lambda_{0j}^{h},\underline{X})$

and obtain the value of $\beta_{j}^{(h+1)}$.

  4. Step 4:

Repeat the steps 2 and 3 for $h=1,2,...$, until converge criteria in Eq. (13)
obtained.

##  3\. Asymptotic Results

The asymptotic properties of the proposed estimators can be derived using
results from empirical process theory. Zhang (2002) proved some results of
about the asymptotic behaviour of the semiparametric pseudo maximum likelihood
estimators when only a single mode of recurrence is observed and later Wellner
and Zhang (2007) modified the results. When recurrence due to multiple modes
are observed, Sreedevi and Sankaran (2020) studied about the asymptotic
properties of cause specific mean functions. We extend the results from
Wellner and Zhang (2000) into a multiple cause scenario and generalize the
results discussed in Sreedevi and Sankaran (2020) to incorporate covariates.
We establish the asymptotic normality and strong consistency of the proposed
estimators.

As we discuss, we estimate $\beta_{j}$ and $\Lambda_{0j}$ for $j=1,2,...,k$ as
the maximum points of the pseudo likelihood function given in Eq.(4). We
assume that the estimators as well as the true value of the parameters include
in the parameter domain $\mathcal{R}\times\mathcal{F}_{j}$ where
$\mathcal{R}\in\mathcal{R}^{d}$ is a bounded convex set and $\mathcal{F}_{j}$
be the class of functions defined as

$\mathcal{F}_{j}\equiv\\{\Lambda_{j}(.):[0,\infty)\to[0,\infty)|~{}~{}\Lambda_{j}(.)~{}~{}\text{is
monotone non-decreasing
with}~{}~{}\Lambda_{j}(0)=0\\}~{}~{}~{}~{}~{}j=1,2,...,k.$

To prove the asymptotic properties of the estimators, we define the following.
Let $\mathcal{B}_{d}$ and $\mathcal{B}$ denote the collection of Borel sets in
$\mathcal{R}_{d}$ and $\mathcal{R}$ respectively. Let $H(.)$ be the
distribution of the covariate vector $Z$ and $\tau$=max$(t)$ and define
$\mathcal{B}_{1}[0,\tau]=\\{B\cap[0,\tau]:B\in\mathcal{B}\\}$ Now we define
the measures $\psi_{j}$, $\eta_{j}$, and $\theta_{j}$ as follows. For
$B,B_{1}\in\mathcal{B}_{1}[0,\tau]$ and $C\in\mathcal{B}_{d}$ define

$\eta_{j}(B\times
C)=\int_{C}\sum\limits_{m=1}^{\infty}P(M=m;J=j|Z=z)\times\sum\limits_{p=1}^{m}P(T_{m,p}\in
B|M=m,Z=z)dH(z).$

A similar measure is defined by Shick and Yu (2000) to study the consistency
of the likelihood estimators for mixed case interval censored data. Define the
$L_{2}$ metric $d_{1}(.)$ in parameter space
$\mathcal{R}\times\mathcal{F}_{j}$ as

$d_{1}((\beta_{j1},\Lambda_{0j1}),(\beta_{j2},\Lambda_{0j2}))=\\{|\beta_{j1}-\beta_{j2}|^{2}+\|\Lambda_{0j1}-\Lambda_{0j2}\|^{2}_{L_{2}(\psi_{1j})}\\}^{\frac{1}{2}},$

where $(\beta_{j1},\Lambda_{0j1})$ and $(\beta_{j2},\Lambda_{0j2})$ are
elements of the parameter space $R^{d}\times\Omega^{+}$ and
$\psi_{j}(B)=\eta_{j}(B\times\mathcal{R}^{d})$.Now to establish the strong
consistency of the estimators, we state the following regularity conditions.
C1: The true parameter values of $\beta_{j}$ and $\Lambda_{0j}$ include in
$\mathcal{R}^{0}\times\mathcal{F}$, where $\mathcal{R}^{0}$ is the interior of
$\mathcal{R}$.
C2: The observation times $T_{M,p}$ are the random variables included in the
bounded interval $[0,\tau]$ for some $\tau\in(0,\infty)$ for all
$p=1,2,...,M$, $M=1,2,...$. Also the measure $\psi_{j}\times H$ on
$([0,\tau]\times\mathcal{R}^{d},\mathcal{B}_{1}[0,\tau]\times\mathcal{B}_{d})$
is absolutely continuous with respect to $\eta_{j}$ for $j=1,2,...,k$ and
$E(M)<\infty$.
C3: For each the true baseline cumulative mean function
$\Lambda_{0j}$,$j=1,2,...,k$ , there exist and $I_{j}\in(0,\infty)$ such that
$\Lambda_{0j}(\tau)\leq I_{j}$.
C4: The function $I_{0j}$ defined as
$I_{0j}(X)\equiv\sum\limits_{p=1}^{M}N_{M,p}\text{log}(N_{M,p})$, satisfies
$P(I_{0j}(X))<\infty$.
C5: The support of $H$, the distribution of covariate vector $Z$ is a bounded
set in $\mathcal{R}^{d}$.
C6: For all $a\in\mathcal{R}^{d},a\neq 0$ and $c\in\mathcal{R}$,
$P(a^{{}^{\prime}}Z\neq c)>0$.
Theorem 1
Under the above regularity conditions C1-C6, and the proposed model specified
by Eq. (1), for every $b<\tau$ such that $\psi_{1j}([b,\tau])>0$,

$d_{1}((\hat{\beta}_{j},\hat{\Lambda}_{0j}I_{[0,b]}),(\beta_{j},\Lambda_{0j}I_{[0,b]}))\to
0~{}~{}~{}~{}~{}~{}~{}\text{a.s}~{}~{}~{}~{}~{}n\to\infty.$

Specifically , when $\psi_{j}({\tau})>0$, we have

$d_{1}((\hat{\beta}_{j},\hat{\Lambda}_{0j}),(\beta_{j},\Lambda_{0j}))\to
0~{}~{}~{}~{}~{}~{}~{}\text{a.s}~{}~{}~{}~{}~{}n\to\infty.$

The proof of the result can be derived by extending the results in Wellner and
Zhang (2007).
To derive the rate of convergence, apart from the above stated regularity
conditions, we also suppose that
C7: For some interval $O_{j}[T]=[\sigma_{j},\tau]$ with some $\sigma_{j}>0$
with $\Lambda_{0j}(\sigma_{j})>0$ and
$P({\bigcap}_{p=1}^{M}){T_{M,p}\in[\sigma_{j},\tau]})=1$.
C8: The number of observations are bounded ie. $P(M\leq m_{0})=1$ for some
$m_{0}\leq\infty$.
C9: For some $z_{0j}\in(0,\infty)$ the function $Z\to
E(e^{z_{0j}N_{j}(\tau)})$ is uniformly bounded for any $Z$ and for all
$j=1,2,..,k$.
C10: There exists a constant $s_{0}>0$ such that $P(T_{M,p}-T_{M,p-1}\geq
s_{0}\text{for all}~{}~{}p=1,2,...,M$)=1. Also, $\psi_{j}(t)$ is absolutely
continuous with respect to a Lebesgue measure with
$0<c_{0}<\psi_{j}^{{}^{\prime}}(t)$ where $c_{0j}$ is a positive constant and
$\psi_{j}^{{}^{\prime}}$ is the derivative of $\psi_{j}$.
C11 : The true baseline cumulative mean functions $\Lambda_{0j}$’s are
differentiable and the derivatives has positive and finite lower and upper
bounds in the observation interval for all $j=1,2,..,k$. ie for each $j$ there
exists a constant $l_{j}$ such that
$\frac{1}{l_{j}}\leq\Lambda_{0j}^{{}^{\prime}}(t)\leq l_{j}\leq\infty$ for
$t\in O[T]$.
Theorem 2:
Under the above stated regularity conditions C7-C11 and the conditions C1-C6
stated to prove consistency of the estimators, for the constant $z_{0j}$
defined in C9, satisfying $z_{0j}\geq 4m_{0}(1+g_{0j})^{2}$ where
$g_{0j}=\sqrt{sign{c_{0j}\Lambda_{0j}^{3}(\sigma_{j})/(24.8l_{j})}}$ and
$\psi_{j}({\tau})>0$,

$n^{\frac{1}{3}}d_{1}(\hat{\beta}_{j},\hat{\Lambda}_{0j}),(\beta_{j},\Lambda_{0j})=O_{p}(1).$

We can see that the rate of convergence of estimators is of order
$n^{-\frac{1}{3}}$ only. Even though the over all convergence rate is
$n^{-\frac{1}{3}}$, we can establish the asymptotic normality of regression
parameters, with the rate of convergence $n^{-\frac{1}{2}}$. Huang (1996)
considered this similar situation for current status data and Sreedevi et al.
(2017) proved similar results for current status data with competing risks.
Theorem 3:
Under the regularity conditions for Theorem 2, the estimator $\hat{\beta}_{j}$
is asymptotically normal and

$n^{\frac{1}{2}}(\hat{\beta}_{j}-\beta)\to_{d}\bar{Z},$

where $\bar{Z}\sim N_{d}(0,\Sigma^{-1}\Theta(\Sigma^{-1})^{{}^{\prime}})$ with

$\Theta=E\Big{(}\sum\limits_{p,p^{{}^{\prime}}=1}^{M}C^{j}_{p,p^{{}^{\prime}}}(Z)[Z-R(M,T_{M,p})][Z-R(M,T_{M,p})]^{{}^{\prime}}\Big{)},$
$\Sigma=E\Big{(}\sum\limits_{p=1}^{M}\Lambda_{j}(T_{M,p})exp(\beta_{j}^{{}^{\prime}}Z)[Z-R(M,T_{M,p})]^{\bigotimes
2}\Big{)},$

in which, $R(M,T_{M,p})\equiv
E(Z\text{exp}(\beta_{j}^{{}^{\prime}}Z)|M,T_{M,p})/E(\text{exp}(\beta_{j}^{{}^{\prime}}Z)|M,T_{M,p})$
and
$C^{j}_{p,p^{{}^{\prime}}}(Z)=\text{Cov}(N_{j}(T_{M,p}),N_{j}(T{M,p^{{}^{\prime}}}))$.
We can see that, in general $\hat{\beta}_{j}$ is not asymptotically efficient,
but when the counts $\\{N_{j}(T_{M,p}),p=1,2,...,M\\}$ consist a cluster of
Poisson count data where the counts within the cluster are independent, the
estimator $\beta_{j}$ become asymptotically efficient for $j=,1,...,k$. Proof
of Theorem 3 can be obtained as a generalisation of the results in Wellner and
Zhang (2007).

##  4\. Simulation Study

We carry out a Monte Carlo simulation study to assess the performance of the
proposed estimation procedure in finite samples. We consider the situation
with two competing risks. The real life situations in reliability and survival
studies are taken as a model to generate panel count data of the form
$\\{M_{i},T_{M_{i},p},N_{M_{i},p}^{1},N_{M_{i},p}^{2},Z_{i}\\}$ for
$p=1,2,...,M_{i}$, $i=1,2,...,n$. We consider
$Z_{i}=\\{Z_{i1},Z_{i2}\\}^{{}^{\prime}}$, as the covariate vector with two
mutually independent components. For each subject, $Z_{i1}$ is generated from
a Bernoulli distribution with probability of success $0.5$ and $Z_{i2}$ is
generated from a Normal distribution with mean $0$ and standard deviation
$0.5$ . The number of observation times $M_{i}$ for each individual is
generated from a discrete uniform distribution $U(1,5)$ for $i=1,2,...,n$.
Thus the maximum number of observations for each individual is restricted upto
5. Then we generated gap times between each observation from uniform
distribution $U(1,5)$. The discrete observation time points $T_{M_{i},p}$ for
$p=1,2,...,M_{i}$ and $i=1,2,...,n$ are generated using the above mentioned
time gaps. Once the observation times are generated, number of recurrences
$\\{N_{M_{i},p}^{1},N_{M_{i},p}^{2}\\}$ are generated from a bivariate Poisson
process given by

$(\Delta N_{M_{i},p}^{1},\Delta
N_{M_{i},p}^{2})\sim\text{BivPo}(\Lambda_{01}(\Delta
T_{M_{i},p})\text{exp}(\beta_{1}^{{}^{\prime}}Z_{i}),\Lambda_{02}(\Delta
T_{M_{i},p})\text{exp}(\beta_{2}^{{}^{\prime}}Z_{i}),\rho),$ (14)

where $\Delta N_{M_{i},p}^{j}=N_{M_{i},p}^{j}-N_{M_{i},p-1}^{j}$ for $j=1,2$;
$\Delta T_{M_{i},p}=T_{M_{i},p}-T_{M_{i},p-1}$, $\Lambda_{01}(t)$ and
$\Lambda_{02}(t)$ are the true baseline functions due to mode l and mode 2,
$\beta_{1}$ and $\beta_{2}$ are the values of regression parameters due to
mode l and mode 2, and $\rho$ is the covariance between the number of
recurrences due to mode 1 and mode 2.

We consider two different forms of $\Lambda_{01}(t)$ and $\Lambda_{02}(t)$ ,
$t$ and $2t$ to generate panel count data. The sample size $n$ takes three
different values $n=50,100,200$. The process is repeated 10000 times to
estimate the efficiency of the estimators. The absolute bias and mean square
error (MSE) of the estimates of $\beta_{1}=\\{\beta_{11},\beta_{12}\\}$ and
$\beta_{2}=\\{\beta_{21},\beta_{22}\\}$ are obtained.

Various parameter values of $\beta_{1}$ and $\beta_{2}$ are considered. Since
the results are similar we present the same only for three different
combinations of $\beta_{1}$ for $\beta_{2}$ in Tables 1-3. To obtain the
convergence we choose $\epsilon=10^{-5}$. The covariance $\rho$ is set to be
$0.5$ in our studies. The simulations are carried out using R programming
language.

Table 1.  Absolute bias and MSE of the estimators of regression coefficients |  | $(\beta_{11},\beta_{12})=(0.5,1)$ | $(\beta_{21},\beta_{22})=(-1,0.5)$
---|---|---|---
True Baseline Function | $n$ | Bias11 | Bias12 | MSE11 | MSE12 | Bias21 | Bias22 | MSE21 | MSE22
$\Lambda_{01}(t)=t,\Lambda_{02}(t)=2t$ | 50 | 0.0248 | 0.0382 | 0.0218 | 0.0108 | 0.1210 | 0.0216 | 0.1080 | 0.0343
100 | 0.0127 | 0.0137 | 0.0099 | 0.0093 | 0.0928 | 0.0199 | 0.0856 | 0.0243
200 | 0.0098 | 0.0117 | 0.0012 | 0.0065 | 0.0720 | 0.0076 | 0.0098 | 0.0105
$\Lambda_{01}(t)=2t,\Lambda_{02}(t)=2t$ | 50 | 0.0198 | 0.0454 | 0.0194 | 0.0121 | 0.0211 | 0.0278 | 0.1097 | 0.0218
100 | 0.0114 | 0.0218 | 0.0089 | 0.0074 | 0.0141 | 0.0124 | 0.0954 | 0.0122
200 | 0.0073 | 0.0098 | 0.0065 | 0.0059 | 0.0069 | 0.0072 | 0.0088 | 0.0076
$\Lambda_{01}(t)=t,\Lambda_{02}(t)=t$ | 50 | 0.0132 | 0.0245 | 0.0110 | 0.0279 | 0.0510 | 0.0199 | 0.0350 | 0.0214
100 | 0.0093 | 0.0135 | 0.0065 | 0.0138 | 0.0061 | 0.0131 | 0.0102 | 0.0166
200 | 0.0043 | 0.0089 | 0.0027 | 0.0020 | 0.0045 | 0.0091 | 0.0071 | 0.0065
Table 2.  Absolute bias and MSE of the estimators of regression coefficients |  | $(\beta_{11},\beta_{12})=(1,0.5)$ | $(\beta_{21},\beta_{22})=(0.5,1)$
---|---|---|---
True Baseline Function | $n$ | Bias11 | Bias12 | MSE11 | MSE12 | Bias21 | Bias22 | MSE21 | MSE22
$\Lambda_{01}(t)=t,\Lambda_{02}(t)=2t$ | 50 | 0.1311 | 0.0393 | 0.0198 | 0.0187 | 0.0312 | 0.0247 | 0.1201 | 0.0298
100 | 0.1008 | 0.0137 | 0.0108 | 0.0082 | 0.0219 | 0.0187 | 0.0916 | 0.0117
200 | 0.0832 | 0.0121 | 0.0084 | 0.0051 | 0.0119 | 0.0931 | 0.0411 | 0.0095
$\Lambda_{01}(t)=2t,\Lambda_{02}(t)=2t$ | 50 | 0.0278 | 0.0298 | 0.0171 | 0.0186 | 0.0217 | 0.0354 | 0.0521 | 0.0221
100 | 0.0211 | 0.0218 | 0.0092 | 0.0110 | 0.0156 | 0.0219 | 0.0954 | 0.0131
200 | 0.0102 | 0.0102 | 0.0059 | 0.0072 | 0.0091 | 0.0141 | 0.0127 | 0.0091
$\Lambda_{01}(t)=t,\Lambda_{02}(t)=t$ | 50 | 0.0538 | 0.0213 | 0.0194 | 0.0232 | 0.0144 | 0.0213 | 0.0212 | 0.0273
100 | 0.0391 | 0.0104 | 0.0102 | 0.0117 | 0.0081 | 0.0121 | 0.0126 | 0.0139
200 | 0.0708 | 0.0091 | 0.0083 | 0.0054 | 0.0042 | 0.0072 | 0.0087 | 0.0076
Table 3.  Absolute bias and MSE of the estimators of regression coefficients |  | $(\beta_{11},\beta_{12})=(1,-2)$ | $(\beta_{21},\beta_{22})=(-1,2)$
---|---|---|---
True Baseline Function | $n$ | Bias11 | Bias12 | MSE11 | MSE12 | Bias21 | Bias22 | MSE21 | MSE22
$\Lambda_{01}(t)=t,\Lambda_{02}(t)=2t$ | 50 | 0.0384 | 0.0421 | 0.0212 | 0.0186 | 0.0492 | 0.0291 | 0.0293 | 0.0418
100 | 0.0276 | 0.0187 | 0.0109 | 0.0115 | 0.0321 | 0.0182 | 0.0172 | 0.0329
200 | 0.0119 | 0.0113 | 0.0083 | 0.0085 | 0.0238 | 0.0732 | 0.0093 | 0.0228
$\Lambda_{01}(t)=2t,\Lambda_{02}(t)=2t$ | 50 | 0.0275 | 0.0471 | 0.0256 | 0.0219 | 0.0221 | 0.0269 | 0.0421 | 0.0253
100 | 0.0262 | 0.0218 | 0.0192 | 0.0143 | 0.0162 | 0.0182 | 0.0321 | 0.0187
200 | 0.0133 | 0.0128 | 0.0102 | 0.0092 | 0.0101 | 0.0116 | 0.0192 | 0.0092
$\Lambda_{01}(t)=t,\Lambda_{02}(t)=t$ | 50 | 0.0218 | 0.0291 | 0.0172 | 0.0267 | 0.0279 | 0.0199 | 0.0401 | 0.0271
100 | 0.0113 | 0.0173 | 0.0108 | 0.0129 | 0.0162 | 0.0131 | 0.0284 | 0.0192
200 | 0.0095 | 0.0121 | 0.0072 | 0.0081 | 0.0119 | 0.0091 | 0.0172 | 0.0125

From simulation studies, we observe that the absolute bias and MSE of the
estimators of regression coefficients approaches zero as sample size
increases. This ensure that the proposed estimators are unbiased with nominal
variance. The minimum value of bias and MSE are obtained when the true base
line function takes the form $\Lambda_{0j}(t)=t$ for both $j=1,2$.

##  5\. Data Analysis

The proposed estimation procedure is applied to a real data on skin cancer
chemoprevention trial given in Sun and Zhao (2013) for illustration. The
primary objective of this study was to evaluate the effectiveness of the drug
DFMO (difluoromethylornithine) in reducing new skin cancers in a population
with a history of non-melanoma skin cancers, basal cell carcinoma and squamous
cell carcinoma. The patients were randomly assigned into two groups a
treatment group with oral DFMO at a daily dose 0.5 gm and a palcebo group with
a matching dosage. The data consist of the details of 290 patients with
history of non-melanoma skin cancers who were supposed to be assessed or
observed every 6 months. However, the real observation and follow up times
differ from patient to patient. The data include the number of recurrences of
two types of recurrent events, basal cell carcinoma (BC) and squamous cell
carcinoma (SC). We treat these two types of cancers as two modes of recurrence
following Sreedevi and Sankaran (2020).

In the data set, the number of observations on an individual varies from 1 to
17 and time of observation varies from 12 to 1766 days. For each individual,
the information on age, gender, DFMO status, number of prior tumours are
observed. We consider all 290 patients in our analysis which include 174 are
males and 116 female. To obtain more explicit conclusions, we analyse the data
on males and females separately by taking the covariate information on DFMO
status and number of prior tumours. Out of 290 patients ,147 were assigned to
the placebo group and reaming 147 were treated with oral DFMO. The number of
prior tumours varies from 1 to 35. The estimates of regression parameters with
corresponding standard errors for males are given in Table 4.

Table 4.  Estimates of the regression parameters with corresponding standard error for males Cause | Covariate | Coefficient | SE | P-value
---|---|---|---|---
BC | DFMO | -0.3715 | 0.2331 | 0.0111
| Prior tumours | 0.0685 | 0.0103 | 0.0005
SC | DFMO | -0.2408 | 0.0460 | 0.0600
| Prior tumours | 0.1013 | 0.0292 | 0.0005

The baseline cause specific cumulative mean functions for males are plotted in
Figure 1. The solid line represents the baseline cumulative mean function for
patients with BC and dotted line represents the baseline cumulative mean
function for patients with SC in Figure 1 and 2.

Figure 1.  Baseline cause specific cumulative mean functions for males

The estimates of regression parameters for females are given in Table 5 and
the baseline cause specific cumulative mean functions are plotted in Figure 2.

Table 5.  Estimates of the regression parameters with corresponding standard error for females Cause | Covariate | Coefficient | SE | P-value
---|---|---|---|---
BC | DFMO | -0.1671 | 0.0347 | 0.0632
| Prior tumours | 0.0666 | 0.0456 | 0.0033
SC | DFMO | 0.9557 | 0.0964 | 0.0320
| Prior tumours | 0.1053 | 0.0458 | 0.0023
Figure 2.  Baseline cause specific cumulative mean functions for females

From Tables 4 and 5 we can see that modes of cancer recurrences basal cell
carcinoma (BC) and squamous cell carcinoma (SC) affect males and females in
different ways. The regression estimators for number of prior tumours is
grater than zero for both males and females and for both modes BC and SC. This
implies that as the number of prior tumours increases hazard rate increases
always. Since the hazard ratio is less than unity, we can say that the drug
DFMO decreases the hazard rate for males with both BC and SC and for females
with BC.

From the plots of baseline cumulative mean functions, we can see that
recurrence rate of BC is higher in males than the recurrence rate of SC upto
1800 days (approximately) and from that point recurrence rate of SC crosses
that of BC, while for females recurrence rate of SC is always lower than that
of BC. The plots also show the difference in recurrence patterns of the events
due BC and SC for males and females.

##  6\. Conclusion

Panel count data with multiple modes of failure often arise in periodic follow
up studies that consider recurrent events exposed to multiple modes. In this
article, we proposed a new proportional mean model for the analysis of panel
count data with multiple modes of recurrence. Estimators for regression
parameters and baseline cumulative mean functions due to each recurrence mode
are derived. A simple iterative procedure is developed for the estimation of
parameters. The finite sample performance of the estimators in terms of bias
and MSE is assessed through a Monte Carlo simulation study. A real data set on
skin cancer chemo prevention trial is analysed using the proposed procedures.

The estimation procedure we developed in this article considered the pseudo
likelihood function of panel count data. Maximum likelihood estimators in this
situation can be developed by extending the results in Wellner and Zhang
(2007), which involves a more complex iterative procedure. An approach based
on estimating equations can also be examined for the regression analysis of
panel count data with multiple recurrence modes. In many situations, rate
functions of the underlying recurrent event process are of importance than
mean functions. Cause specific rate functions developed by Sankaran et al.
(2021) can be used to study panel count data, when subjects are exposed to
multiple recurrence modes.

## References

  * [1] Chiou, S. H., Xu, G., Yan, J. and Huang, C. Y. (2018). Semiparametric estimation of the accelerated mean model with panel count data under informative examination times, Biometrics, 74(3), 944–953.
  * [2] Chiou, S. H., Huang, C. Y., Xu, G. and Yan, J. (2019). Semiparametric regression analysis of panel count data: A practical review, International Statistical Review, 87(1), 24–43.
  * [3] Huang, J. (1996). Efficient estimation for the proportional hazards model with interval censoring, The Annals of Statistics, 24(2), 540–568
  * [4] Huang, C. Y., Wang, M. C. and Zhang, Y. (2006). Analysing panel count data with informative observation times, Biometrika, 93(4), 763–775.
  * [5] Hu, X. J., Sun, J. and Wei, L. J. (2003). Regression parameter estimation from panel counts, Scandinavian Journal of Statistics, 30(1), 25–43.
  * [6] Kalbfleisch, J. D. and Lawless, J. F. (1985). The analysis of panel data under a Markov assumption, Journal of the American Statistical Association, 80(392), 863–871.
  * [7] Kim, Y. J. (2007). Analysis of panel count data with measurement errors in the covariates, Journal of Statistical Computation and Simulation, 77(2), 109–117.
  * [8] Lawless, J. F. and Zhan, M. (1998). Analysis of interval grouped recurrent event data using piecewise constant rate functions, Canadian Journal of Statistics, 26(4), 549–565.
  * [9] Sankaran, P. G., Ashlin Mathew, P. M. and Sreedevi E. P. (2021). Cause specific rate functions for panel count data with multiple modes of recurrence, Journal of Indian Statistical Association, 58 (2), In press.
  * [10] Schick, A. and Yu, Q. (2000). Consistency of the GMLE with mixed case interval censored data, Scandinavian Journal of Statistics, 27(1), 45–55.
  * [11] Sreedevi, E. P., Sankaran, P. G. and Dewan, I. (2017). A semi-parametric regression model for current status competing risks data, Journal of the Indian Statistical Association, 55(1), 35–61.
  * [12] Sreedevi, E. P. and Sankaran, P. G. (2020). Nonparametric inference for panel count data with competing risks, Journal of Applied Statistics, 1–14 DOI:10.1080/02664763.2020.1795816.
  * [13] Sun, J. and Kalbfleisch, J. D. (1995). Estimation of the mean function of point processes based on panel count data, Statistica Sinica, 5(1), 279–290.
  * [14] Sun, J. and Wei, L. J. (2000). Regression analysis of panel count data with covariate-dependent observation and censoring times, Journal of the Royal Statistical Society: Series B (Statistical Methodology), 62(2), 293–302.
  * [15] Sun, J. (2007). The Statistical Analysis of Interval-Censored Failure Time Data. Springer, New York.
  * [16] Sun, J., Tong, X. and He, X. (2007). Regression analysis of panel count data with dependent observation times, Biometrics, 63(4), 1053–1059.
  * [17] Sun, J. (2009) Panel count data. Handbook of Statistical methods in Life and Health Sciences, Editor: Balakrishnan, N., John Wiley and Sons Ltd.
  * [18] Sun, J. and Zhao, X. (2013). Statistical Analysis of Panel Count Data. Springer, New York.
  * [19] Thall, P. F. (1988). Mixed Poisson likelihood regression models for longitudinal interval count data, Biometrics, 44(1), 197–209.
  * [20] Thall, P. F. and Lachin, J. M. (1988). Analysis of recurrent events: Nonparametric methods for random-interval count data, Journal of the American Statistical Association, 83(402), 339–347.
  * [21] Wellner, J. A. and Zhang, Y. (2000). Two estimators of the mean of a counting process with panel count data, The Annals of Statistics, 28(3), 779–814.
  * [22] Wellner, J. A. and Zhang, Y. (2007). Two likelihood-based semiparametric estimation methods for panel count data with covariates, The Annals of Statistics, 35(5), 2106–2142.
  * [23] Xu, D., Zhao, H. and Sun, J. (2018). Joint analysis of interval censored failure time data and panel count data, Lifetime Data Analysis, 24(1), 94–109.
  * [24] Zhang, Y. (2002). A semiparametric pseudo likelihood estimation method for panel count data, Biometrika, 89(1), 39–48.
  * [25] Zhao, X. and Tong, X. (2011). Semiparametric regression analysis of panel count data with informative observation times, Computational Statistics and Data Analysis, 55(1), 291–300.


††thanks:  This work is supported by Knut and Alice Wallenberg Foundation and
the Swedish Research Council.
$*$ Corresponding author.
E-mail addresses<EMAIL_ADDRESS>(C. Wu<EMAIL_ADDRESS>(M. Zhang<EMAIL_ADDRESS>(D. V. Dimarogonas)

, ,

# On Disturbance Propagation in Vehicular Platoons with Different
Communication Ranges

Chengshuai Wu∗ Meng Zhang  Dimos V. Dimarogonas  School of Electrical
Engineering and Computer Science, KTH Royal Institute of Technology and
Digital Futures, Stockholm SE-100 44, Sweden.  School of Cyber Science and
Engineering, Xi’an Jiaotong University, Xi’an 710049, China.

###### Abstract

In the control of vehicular platoons, the disturbances acting on one vehicle
can propagate and affect other vehicles. If the disturbances do not amplify
along the vehicular string, then it is called string stable. However, it is
usually difficult to achieve string stability with a distributed control
setting, especially when a constant spacing policy is considered. This note
considers the string unstable cases and studies disturbance propagation in a
nonlinear vehicular platoon consisting of $n+1$ vehicles where the (virtual)
leading vehicle provides the reference for a constant spacing policy. Apart
from the communications between consecutive vehicles, we also assume that each
vehicle can receive information from $r$ neighbors ahead, that is, the
vehicular platoon has communication range $r$. For the maximal overshoot of
the inter-vehicular spacing errors, we explicitly show that the effect of
disturbances, including the external disturbances acting on each vehicle and
the acceleration of the leading vehicle, is scaled by
$O(\sqrt{\left\lceil\frac{n}{r}\right\rceil})$ for a fixed $n$. This implies
that disturbance propagation can be reduced by increasing communication range.
Numerical simulation is provided to illustrate the main results.

###### keywords:

Platoon, contraction, matrix measure, input-to-state stability, string
stability.

##  1 Introduction

Vehicular platooning is an intelligent transportation application which allows
a string of vehicles to move in a closely spaced manner. This can decrease the
aerodynamic drag and increase the traffic throughput (Varaiya, 1993).
Centralized control is impractical when the vehicular platoons are large in
length considering the limited wireless communication range and the delay in
transmission. Therefore, distributed control protocols are considered (see
e.g., Zheng et al. (2015); Herman et al. (2017); Peters et al. (2014)).
However, the distributed design leads to the so called disturbance propagation
in vehicular platoons.

Many works consider a safety and performance criteria related to disturbance
propagation called _string stability_ (see the recent see the recent tutorial
(Feng et al., 2019)). There are various definitions for string stability in
existing works, but they all assert similar properties that the effect of
disturbances does not amplify along the vehicular string and is independent of
the platoon length. Depending on where the disturbances under concern act on,
string stability w.r.t. the external input (e.g., acceleration) of the leading
vehicle (Ploeg et al., 2013), and string stability w.r.t. the external
disturbances on each vehicles (Besselink and Johansson, 2017; Monteil et al.,
2019) can be defined. Depending on which kind of signal norms is utilized to
characterize the effect of the disturbances, there are mainly two notions,
i.e., $\mathcal{L}_{2}$ string stability (Ploeg et al., 2013; Stüdli et al.,
2018; Herman et al., 2016) and $\mathcal{L}_{\infty}$ string stability
(Swaroop and Hedrick, 1996; Besselink and Johansson, 2017; Monteil et al.,
2019). Specifically, $\mathcal{L}_{2}$ string stability can be specified using
$H_{\infty}$ norms of the transfer functions related to the dynamics of each
vehicles. The characterization of $\mathcal{L}_{\infty}$ string stability
usually relies on Lyapunov methods (Swaroop and Hedrick, 1996; Besselink and
Johansson, 2017). Recently, Monteil et al. (2019) shows that
$\mathcal{L}_{\infty}$ string stability can be elegantly specified based on
contraction theory (Aminzare and Sontag, 2014). This method is extended in
(Silva et al., 2021) where a disturbance rejection control design is proposed.

Despite the nice properties ensured by string stability, there are many
negation results suggest that it is usually difficult to achieve string
stability especially when the constant spacing policy is considered. For
example, Seiler et al. (2004); Barooah and Hespanha (2005) show that vehicular
platoons with double integrator dynamics using only relative spacing
information and constant spacing policies are always $\mathcal{L}_{2}$ string
unstable for any linear controller. Yadlapalli et al. (2006) shows that string
stability cannot be guaranteed for non-cyclic vehicular platoons with constant
spacing policies if the communication range is limited. Peters et al. (2016)
shows that cyclic vehicular platoons with constant spacing policies using only
the information of the predecessor are not even stable when the platoon length
is large enough. Recently, Farnam and Sarlette (2019) provide an impossibility
result for string stability which covers a broad range of control designs. It
is known that string stability can be achieved if the leading vehicle can
broadcast its information to all the other vehicles (see e.g., Monteil et al.
(2019); Besselink and Johansson (2017)). However, it is not practical when the
platoon length is large and the broadcast delays should be considered (Peters
et al., 2014). Another solution to string stability is using time headway
policies (see e.g. Ploeg et al. (2013)) which rely on absolute velocity
information.

For the string unstable cases, analyzing the effect of disturbance propagation
is of both theoretical and practical interest. For linear vehicular platoons,
this can be analyzed by a transfer function matrix that relates the vector of
all external disturbances acting on the platoon to the vector of all spacing
errors (see e.g., Seiler et al. (2004); Yadlapalli et al. (2006)). Lyapunov
methods for studying disturbance propagation have also been reported (see
e.g., Li et al. (2019); Herman et al. (2017)), and the results are essentially
on the convergence rate of the platoon to a desired spacing configuration
since disturbance propagation can be reduced by increasing the convergence
rate.

Unlike the aforementioned works, the disturbance propagation problem studied
in this paper is formulated by the relation between the disturbances and the
maximal overshoot of the inter-vehicular spacing errors. Note that if this
relation is independent of the platoon length, then a $\mathcal{L}_{\infty}$
string stability-like result is achieved. Similar to (Yadlapalli et al., 2006;
Li et al., 2019), we study vehicular platoons with variable communication
ranges, that is, each vehicle can receive information from certain amount of
neighbors ahead. Our main contribution is related to a contraction theory
framework which allows to analyze nonlinear vehicular platoons and can provide
explicit bound for the maximal overshoot of spacing errors subject to external
disturbances acting on each vehicle and the acceleration of the leading
vehicle. For a vehicular platoon consisting of $n+1$ vehicles with
communication range $r$, we show that the effect of disturbances on the
aforementioned bound is scaled by
$O(\sqrt{\left\lceil\frac{n}{r}\right\rceil})$ for a fixed $n$.

The remainder of this paper is organized as follows. The next section reviews
some basic notions and results in contraction theory. Section 3 presents the
configuration and modeling of the studied vehicular platoons and a unified
control protocol for different communication ranges. The main result is given
in Section 4. Section 5 describes a byproduct of the main result, i.e., a
$\mathcal{L}_{\infty}$ string stability result. Simulation analysis is
provided in Section 6.

##  2 Preliminaries

Let $|\cdot|:\mathbb{R}^{n}\to\mathbb{R}_{\geq 0}$ denote a vector norm. The
induced matrix norm $||\cdot||:\mathbb{R}^{n\times p}\to\mathbb{R}_{\geq 0}$
is $||A||:=\max_{|x|=1}|Ax|$. The associated matrix measure
$\mu:\mathbb{R}^{n\times n}\to\mathbb{R}$ is

$\mu(A):=\lim_{\varepsilon\downarrow 0}\frac{||I_{n}+\varepsilon
A||-1}{\varepsilon}.$

For $p\in\\{1,2,\infty\\}$, denote $\mu_{p}$ as the matrix measure associated
with the $p$-norm $|\cdot|_{p}$. For any measurable function
$\theta:\mathbb{R}_{\geq 0}\to\mathbb{R}^{n}$ and a vector norm $|\cdot|$,
define $||\theta||:=\operatorname{ess.sup}_{t\geq 0}|\theta(t)|$, i.e., the
signal $\mathcal{L}_{\infty}$ norm induced by the vector norm $|\cdot|$. For
example, we use $||\theta||_{\infty}$ to denote the signal
$\mathcal{L}_{\infty}$ norm induced by $|\cdot|_{\infty}$. For a signal
$u(t)\in\mathbb{R}^{n}$ defined on $[0,T)$ and for some $\tau\in[0,T)$, define
the related truncation signals:

$\displaystyle u_{\tau}(t)=\begin{cases}u(t),&t\in[0,\tau],\\\
0,&t\in(0,\infty).\end{cases}$ $\displaystyle
u^{\tau}(t)=\begin{cases}0,&t\in[0,\tau),\\\ u(t),&t\in[\tau,T).\end{cases}$

Then, for $0<\tau_{1}<\tau_{2}<T$, a truncation on $t\in[\tau_{1},\tau_{2}]$
of $u(t)$ is denoted $u_{\tau_{2}}^{\tau_{1}}$.

Next we review some basic results on contractive systems (Aminzare and Sontag,
2014). Consider the nonlinear time-varying system:

$\dot{x}(t)=f(t,x(t)),$ (1)

where $f:\mathbb{R}_{+}\times\Omega\to\mathbb{R}^{n}$ is $C^{1}$ and $\Omega$
is a convex subset of $\mathbb{R}^{n}$. For $t\geq t_{0}\geq 0$ and
$x_{0}\in\Omega$, we assume that (1) admits a unique solution
$x(t,t_{0},x_{0})$ for all $t\geq t_{0}$, and that $x(t,t_{0},x_{0})\in\Omega$
for all $t\geq t_{0}$. From here on we always take $t_{0}=0$ and write
$x(t,x_{0})$ for $x(t,0,x_{0})$.

The system (1) is called contractive if there exist a vector norm
$|\cdot|:\mathbb{R}^{n}\to\mathbb{R}_{\geq 0}$ and $\eta>0$ such that for any
$a,b\in\Omega$,

$|x(t,a)-x(t,b)|\leq\exp(-\eta t)|a-b|,\text{ for all }t\geq 0.$ (2)

In other words, any two trajectories contract to each other at an exponential
rate. Hence, all the trajectories eventually converge to a unique steady-state
solution.

Let $J(t,x):=\frac{\partial}{\partial x}f(t,x)$ denote the Jacobian of the
vector field $f$ with respect to $x$. A sufficient condition (Russo et al.,
2010) guaranteeing (2) is that the matrix measure
$\mu(\cdot):\mathbb{R}^{n\times n}\to\mathbb{R}$, associated with $|\cdot|$,
satisfies

$\mu(J(t,x))\leq-\eta,\text{ for all }t\geq 0,x\in\Omega.$ (3)

The next result given by Desoer and Haneda (1972) shows that a contractive
system with an additive input is input-to-state stable (ISS) (Sontag et al.,
1989, Def. 2.1).

######  Proposition 1.

Consider the system (1) with an input $u(t)\in\mathbb{R}^{n}$, i.e.,

$\dot{x}(t)=f(t,x(t))+u(t),$ (4)

where $u:\mathbb{R}_{+}\to\mathbb{R}^{n}$ is piecewise continuous. Assume that
condition (3) holds with $\Omega=\mathbb{R}^{n}$, and $f(t,0)=0$ for all
$t\geq 0$. Then,

$|x(t,x_{0})|\leq\exp(-\eta
t)|x_{0}|+\int_{0}^{t}\exp(-\eta(t-s))|u(s)|\mathop{}\\!\mathrm{d}s,$ (5)

for all $t\geq 0$. That is, there exists a class $\mathcal{KL}$ function
$\beta:\mathbb{R}_{\geq 0}\times\mathbb{R}_{\geq 0}\to\mathbb{R}_{\geq 0}$ and
a class $\mathcal{K}$ function $\gamma:\mathbb{R}_{\geq 0}\to\mathbb{R}_{\geq
0}$ such that

$|x(t,x_{0})|\leq\beta(|x_{0}|,t)+\gamma(||u_{t}||),\text{ for all }t\geq 0,$
(6)

where $\beta(s,t)=\exp(-\eta t)s$ and $\gamma(s)=s/\eta$.

##  3 Problem formulation

Consider the longitudinal formation control problem of a vehicular string
consisting of $(n+1)$ vehicles. The position and speed of the $i$th vehicle is
denoted $p_{i}(t)$ and $v_{i}(t)$, $i=0,1,\dots,n$, respectively. The state of
leading vehicle $(p_{0}(t),v_{0}(t))$ provides a reference for the platoon.
The dynamics of the $i$th vehicle, $i=1,\dots,n$, is described as a second
order linear system:

$\displaystyle\dot{p}_{i}$ $\displaystyle=v_{i},$ (7) $\displaystyle
m_{i}\dot{v}_{i}$ $\displaystyle=u_{i}+\theta_{i},$

where $m_{i}$ is the mass of the $i$th vehicle, $\theta_{i}(t)$ denote a
disturbance, and $u_{i}$ is the control input of the $i$th vehicle. Note that
the longitudinal dynamics of vehicles is typically modeled as a nonlinear
system, and (7) can be obtained via feedback linearization (see e.g. Stankovic
et al. (2000)).

The control goal is to maintain certain designed inter-vehicular distances for
the platoon, i.e., the constant spacing policy. To formalize this, define a
new state $x\in\mathbb{R}^{n}$ whose $i$th entry is

$x_{i}:=p_{i-1}-p_{i},\quad i=1,\dots,n,$

i.e., the $i$th inter-vehicular distance. Let $e\in\mathbb{R}^{n}$ be a vector
with positive entries, and its $i$th entry $e_{i}$, $i=1,\dots,n$, denote the
desired $i$th inter-vehicular distance. We will exploit a set of decentralized
tracking controllers such that the $i$th vehicle can track a desired velocity
given by

$v_{di}(t):=d_{i}(x_{i},x_{i+1})+v_{0},~{}i=1,\dots,n.$ (8)

We set $x_{i+1}=0$, i.e., $d_{n}(x_{n},x_{n+1})=d_{n}(x_{n})$. The mappings
$d_{i}$, $i=1,\dots,n$, are $C^{1}$ w.r.t the arguments and satisfies
$d_{i}(e_{i},e_{i+1})=0$. That is, the velocity of the platoon is $v_{0}(t)$
when the desired configuration $x=e$ is achieved. Throughout the note we
consider a set of $d_{i}$’s such that there exist positive constants
$\eta_{1}$ and $c$,

$\displaystyle\frac{\partial d_{i}}{\partial x_{i}}>0,~{}\frac{\partial
d_{i}}{\partial x_{i+1}}<0,~{}\frac{\partial d_{i}}{\partial
x_{i}}+\frac{\partial d_{i}}{\partial x_{i+1}}\geq\eta_{1},$ (9)
$\displaystyle\max\left\\{\left|\frac{\partial d_{i}}{\partial
x_{i}}\right|,\left|\frac{\partial d_{i}}{\partial
x_{i+1}}\right|\right\\}\leq c,$ (10)

for all $i=1,\dots,n$, and all $x\in\mathbb{R}^{n}$.

######  Remark 1 (Design of $d_{i}$’s).

The mappings $d_{i}(x_{i},x_{i+1})$ can be viewed as a formation protocol
(Mesbahi and Egerstedt, 2010) for the platoon with the first order dynamics
$p_{i}=v_{i}$, $i=1,\dots,n$. For example, we can design $d_{i}$’s as the
standard linear formation protocol

$d_{i}(x_{i},x_{i+1})=\ell_{i}^{p}(x_{i}-e_{i})-\ell_{i}^{f}(x_{i+1}-e_{i+1}),$
(11)

where $\ell^{p}_{i},\ell_{i}^{f}>0$, $i=1,\dots,n$. Note that Conditions
(9)-(10) hold in this case with $\ell_{i}^{p}>\ell_{i}^{p}$ for all $i$.

We assume that consecutive vehicles can communicate information including
relative positions and velocities. Additionally, each vehicle can receive
information from $r$ neighboring vehicles ahead. In other words, the platoon
has communication range $r\in\\{1,\dots,n\\}$.

We propose the following control protocol:

$\begin{split}\frac{1}{m_{i}}u_{i}&=-k_{i}\left(v_{i}-\sum\limits_{j=0}^{r-1}d_{i-j}-v_{i-r}\right)\\\
&~{}+\frac{\partial}{\partial x_{i}}d_{i}(x_{i},x_{i+1})(v_{i-1}-v_{i})\\\
&~{}+\frac{\partial}{\partial
x_{i+1}}d_{i}(x_{i},x_{i+1})(v_{i}-v_{i+1}),~{}i=1,\dots,n,\end{split}$ (12)

where $k_{i}>0$, $i=1,\dots,n$, are the control gains.

The next section details the main results, which provide an upper bound for
the maximal overshoot of the inter-vehicular spacing errors, i.e.,
$\max_{i}\\{||x_{i}-e_{i}||\\}$.

##  4 Main results

We first show that the closed-loop systems of (7) with the control designs
(12) can reformulated as a singular perturbed systems (see e.g. (Khalil, 2002,
Ch. 11)).

Define a new state $z\in\mathbb{R}^{n}$ with its $i$th entry as

$z_{i}:=v_{i}-v_{di}.$

Let $k_{min}:=\min_{i}\\{k_{i}\\}$ and write the control gains $k_{i}$’s as

$k_{i}=\bar{k}_{i}k_{\operatorname{min}},~{}i=1,\dots,n.$

That is, $\bar{k}_{i}$ can be viewed as a relative control gain w.r.t.
$k_{min}$, and satisfies $\bar{k}_{i}\geq 1$ for all $i\in\\{1,\dots,n\\}$.
Define a parameter $\epsilon:=\frac{1}{k_{min}}$ and let $d_{i}=0$ and
$z_{i}=0$ for all $i\leq 0$. Then the closed-loop system of (7) with the
control laws (12) is

$\begin{split}\dot{x}_{i}&=d_{i-1}-d_{i}+z_{i-1}-z_{i},\\\
\epsilon\dot{z}_{i}&=-\bar{k}_{i}\left(z_{i}-z_{i-r}-\sum\limits_{j=1}^{r}d_{i-j}\right)+\epsilon\left(\frac{\theta_{i}}{m_{i}}-\dot{v}_{0}\right),\\\
\end{split}$ (13)

with $i=1,\dots,n$. For simplicity, we define the overall disturbance signal
$w\in\mathbb{R}^{n}$ whose $i$th entry is

$w_{i}=\frac{\theta_{i}}{m_{i}}-\dot{v}_{0},~{}i=1,\dots,n.$

Note that the parameter $\epsilon$ can be made arbitrarily small by exploiting
a high gain design. This justifies that the above system is indeed a
singularly perturbed system whose stability usually depends on the existence
of a small enough $\epsilon$.

By the standard procedure of singular perturbation, i.e., time scale
decomposition, system (13) can be decomposed into the so-called _reduced
system_ and _boundary layer system_. Specifically, let $\epsilon=0$, then the
$z$-system achieves the steady state $h(x)\in\mathbb{R}^{n}$ instantly defined
by

$h_{i}-h_{i-r}-\sum\limits_{j=1}^{r}d_{i-j}=0,~{}i=1,\dots,n,$ (14)

where $h_{i}(x)$ denote the $i$th entry of $h(x)$. A straightforward
calculation yields

$h_{i}(x)=\sum_{j=1}^{i-1}d_{j}(x_{j},x_{j+1}),~{}i=1,\dots,n.$ (15)

Define the error state

$y:=z-h(x),$

and a fast time scale $\tau:=\frac{t}{\epsilon}$. Then the system (13) in the
$(x,y)$-coordinate is

$\displaystyle\dot{x}$ $\displaystyle=f(x)+Gy,$ (16)
$\displaystyle\frac{\mathop{}\\!\mathrm{d}y}{\mathop{}\\!\mathrm{d}\tau}$
$\displaystyle=Ay-\epsilon\frac{\partial h}{\partial x}(f(x)+Gy)+\epsilon w,$
(17)

where $f:=\begin{bmatrix}-d_{1}&-d_{2}&\cdots&-d_{n}\end{bmatrix}^{T}$. Let
$G_{ij}$ and $A_{ij}$ denote the $ij$th entry of the $n\times n$ matrices $G$
and $A$, respectively. Then,

$\displaystyle G_{ij}$ $\displaystyle=\begin{cases}-1,&i=j,\\\ 1,&i=j+1,\\\
0,&\text{otherwise}.\end{cases}$ (18) $\displaystyle A_{ij}$
$\displaystyle=\begin{cases}-\bar{k}_{i},&i=j,\\\ \bar{k}_{i},&i=j+r,\\\
0,&\text{otherwise}.\end{cases}$ (19)

Before the main theorem, we first give the next two instrumental results
studying the matrix measures of the Jacobians of the two systems (16) and
(17), i.e., contractivity of (16) and (17). The next result shows that the
systems (16) is contractive if $y$ is treated as a parameter.

######  Lemma 1.

Define $J(x)=\frac{\partial}{\partial x}f(x).$ Condition (9) ensures that
$\mu_{\infty}(J(x))\leq-\eta_{1}$ for all $x\in\mathbb{R}^{n}$.

The proof is based on that $\mu_{\infty}(J)=\max_{i}\\{J_{ii}+\sum_{j\neq
i}|J_{ij}|\\}$.

The next result gives a sufficient condition such that the system (17) is
contractive when $x$ is viewed as a parameter. Here, we denote $P\leq Q$ for
$P,Q\in\mathbb{R}^{n\times n}$ if $P_{ij}\leq Q_{ij}$ for all
$i,j\in\\{1,\dots,n\\}$. For a symmetric matrix $P\in\mathbb{R}^{n\times n}$,
we denote $P\succ[\succeq]0$ if $P$ is positive [semi]definite.

######  Lemma 2.

Define $\bar{J}(x)=A-\epsilon\frac{\partial}{\partial x}h(x)G.$ Given
$a\in\mathbb{R}^{n}$ and decompose it as

$a=\begin{bmatrix}a^{1}\\\ \vdots\\\
a^{m}\end{bmatrix},~{}a^{i}\in\mathbb{R}^{r},~{}i=1,\dots,m-1,a^{m}\in\mathbb{R}^{n-(m-1)r},$

where $m:=\lceil\frac{n}{r}\rceil$. Let $D\in\mathbb{R}^{m\times m}$ be a
positive definite diagonal matrix. Define a vector norm
$|\cdot|_{*}:\mathbb{R}^{n}\to\mathbb{R}_{\geq 0}$ by

$|a|_{*}:=\left|\begin{bmatrix}|a^{1}|_{\infty}\\\ \vdots\\\
|a^{m}|_{\infty}\end{bmatrix}\right|_{D,2},$ (20)

where $|\cdot|_{D,2}:\mathbb{R}^{m}\to\mathbb{R}_{\geq 0}$ denotes a
$D$-weighted 2-norm defined by $|b|_{D,2}:=|Db|_{2}$ for $b\in\mathbb{R}^{m}$.
Let $\mu_{*}$ denote the matrix measure induced by $|\cdot|_{*}$. For any
$\bar{\epsilon}>0$, Conditions (9) and (10) ensure that if

$-\bar{k}_{i}+2\bar{\epsilon}c(r-1)<0,~{}i=1,\dots,n,$ (21)

then there exist $\eta_{2}>0$ and a $D$ such that

$\mu_{*}(\bar{J}(x))\leq-\eta_{2},\text{ for all
}\epsilon\in(0,\bar{\epsilon}],\text{ and all }x\in\mathbb{R}^{n}.$ (22)

###### Proof.

The arguments are based a result by Ström (1975) (see also (Russo et al.,
2013, Thm. 2), (Ofir et al., 2022a, Thm. 1)).

Consider a $m\times m$ partition of $\bar{J}(x)$ where the $ij$th block is
denoted $\bar{J}^{ij}(x)\in\mathbb{R}^{n_{i}\times n_{j}}$ with $n_{i}=r$ for
$i\in\\{1,\dots,m-1\\}$, and $n_{m}=n-r(m-1)$. Define
$B(x)\in\mathbb{R}^{m\times m}$, with its $ij$th entry denoted $B_{ij}(x)$, by

$B_{ij}(x):=\begin{cases}\mu_{\infty}(\bar{J}^{ij}(x)),&i=j,\\\
||\bar{J}^{ij}(x)||_{\infty},&i\neq j.\end{cases}$ (23)

Note that $B(x)$ is Metzler. Condition (10) ensures that $B(x)$ is bounded for
any $x\in\mathbb{R}^{n}$. That is, there exists a constant Metzler matrix
$\bar{B}\in\mathbb{R}^{n\times n}$ such that $B(x)\leq\bar{B}$ for all
$x\in\mathbb{R}^{n}$. Let $\mu_{D,2}$ denote the matrix measure induced by
$|\cdot|_{D,2}$. Ström (1975) shows that

$\mu_{*}(\bar{J}(x))\leq\mu_{D,2}(B(x))\leq\mu_{D,2}(\bar{B}),\text{ for all
}x\in\mathbb{R}^{n}.$ (24)

It is worth pointing out that $|\cdot|_{D,2}$ is a monotonic norm, which is
necessary for (24) to hold.

A straightforward calculation shows that the matrix
$F(x):=-\epsilon\frac{\partial}{\partial x}h(x)G$ is a lower triangular matrix
satisfying

$|F_{ij}(x)|\leq 2\epsilon c,\text{ for all }i>j.$ (25)

Recall (19), then the matrices $\bar{J}(x)$ and $B(x)$ are also lower
triangular. The diagonal entries of $\bar{J}(x)$ are

$-\bar{k}_{1},~{}\epsilon\frac{\partial d_{1}}{\partial
x_{2}}-\bar{k}_{2},\dots,~{}\epsilon\frac{\partial d_{n-1}}{\partial
x_{n}}-\bar{k}_{n}.$

Note that they are all negative since for any $\epsilon>0$ due to Condition
(9). Therefore, Condition (21) ensures that there exists $c_{1}>0$ such that
if $\epsilon\in(0,\bar{\epsilon}]$, then

$\mu_{\infty}(\bar{J}^{ii}(x))<-c_{1},\text{ for all
}i\in\\{1,\dots,m\\},~{}x\in\mathbb{R}^{n}.$

Since $B(x)$ is lower triangular, its eigenvalues are the diagonal entries,
i.e., $\mu_{\infty}(\bar{J}^{ii}(x))$, $i=1,\dots,m$. Therefore, $B(x)$ is
Hurwitz, and thus there exists a lower triangular and Hurwitz $\bar{B}$. It is
well known that a Hurwitz and Metzler matrix is diagonally stable, i.e., there
exists a diagonal matrix $D\succ 0$ and $\eta_{2}>0$ such that
$D^{2}\bar{B}+\bar{B}^{T}D^{2}\preceq-2\eta_{2}I_{m}$. This implies that
$\mu_{D,2}(B)=\mu_{2}(DBD^{-1})\leq-\eta_{2}$. Then (22) follows from (24) and
(25). $~{}~{}~{}~{}\square$

Note that the parameter $\eta_{2}$, i.e., the contraction rate of the system
(17), depends on the communication range $r$ implicitly since the vector norm
$|\cdot|_{*}$ in (20) is related to $r$. The next example suggests that a
larger $\eta_{2}$ can be achieved by increasing $r$.

######  Example 1.

Let $k_{1}=\cdots=k_{n}$, i.e., $\bar{k}_{1}=\cdots=\bar{k}_{n}=1$. Consider
the extreme case $\epsilon=0$, then $\bar{J}(x)=A$ and the matrix
$B(x)\in\mathbb{R}^{n\times n}$ in (23) is given by

$B_{ij}:=\begin{cases}-1,&i=j,\\\ 1,&i=j+1,\\\
0,&\text{otherwise}.\end{cases}$

The symmetric part of $B$, i.e., $B_{\operatorname{sym}}:=(B+B^{T})/2$, is an
$m\times m$ tridiagonal matrix with $-1$ on the main diagonal, and $1/2$ on
the sub- and super-diagonals. Yueh (2005) shows that the eigenvalues of such a
matrix are

$\lambda_{i}=-1+\cos\Big{(}\frac{i\pi}{m+1}\Big{)},~{}i=1,\dots,m.$

Therefore,

$\mu_{2}(B)=\lambda_{\operatorname{max}}(B_{\operatorname{sym}})=-1+\cos\left(\frac{\pi}{m+1}\right)<0.$
(26)

That is, the matrix $D$ in this case can be selected as $I_{m}$. Eq. (26)
implies that $\mu_{2}(B)$ is increasing w.r.t. $m$ if $\epsilon=0$. For a
small enough $\epsilon>0$, this property still holds since $B(x)$ depends on
$\epsilon$ continuously and $\mu_{2}(B(x))$ is continuous w.r.t. the entries
of $B(x)$. Recall that $m:=\lceil\frac{n}{r}\rceil$. This shows that a larger
communication range $r$ can lead to a smaller $\mu_{2}(B)$, that is, a larger
$\eta_{2}$.

######  Remark 2.

Note that Lemmas 1 and 2 show that the system (16)-(17) can be viewed as an
interconnection of two contracive systems, or two ISS systems by Proposition
1. It is well known that the small gain theorem can be applied in this case.
However, it requires that the term $\frac{\partial h}{\partial x}f(x)$ in (17)
is bounded for all $x\in\mathbb{R}^{n}$. For the incremental version of small
gain theorem (Russo et al., 2013), it requires that $\frac{\partial}{\partial
x}\left(\frac{\partial h}{\partial x}f(x)\right)$ is bounded for all
$x\in\mathbb{R}^{n}$. Note that these two bounded conditions are not required
in this work.

Based on Lemmas 1 and 2, the main theorem below shows that the platooning
system (7) with the control protocol (12) is semi globally stable. The
argument is related to finding a small enough $\epsilon$ (i.e., a large enough
$k_{\operatorname{min}}$). It is inspired by (Christofides and Teel, 1996,
Thm. 1), which is a Lyapunov method based on implicit constructions. In our
case, these implicit constructions can be made explicit by contraction theory.
To facilitate the analysis, define

$\begin{split}||G||_{\infty,*}:=\max\limits_{a\neq
0}\frac{|Ga|_{\infty}}{|a|_{*}},~{}a\in\mathbb{R}^{n}.\end{split}$

Note that $|a|_{*}\geq|a|_{\infty}$ for all $a\in\mathbb{R}^{n}$. Hence,

$||G||_{\infty,*}\leq\max\limits_{a\neq
0}\frac{|Ga|_{\infty}}{|a|_{\infty}}=|G|_{\infty}=2.$ (27)

######  Theorem 2.

Consider the $(x,y)$-system (16)-(17) where the mappings $d_{i}$’s satisfy
Conditions (9)-(10), and the parameters $\bar{k}_{i}$’s satisfy Condition (21)
for some $\bar{\epsilon}>0$. For each pair of positive constants
$(\delta,\sigma)$, there exists $\epsilon^{*}\in(0,\bar{\epsilon}]$, such that
if

$\max\\{|x(0)-e|_{\infty},|y(0)|_{*},||w||_{*}\\}\leq\delta$

and $\epsilon\in(0,\epsilon^{*}]$, then

$\displaystyle|x(t)-e|_{\infty}$
$\displaystyle\leq\exp(-\eta_{1}t)|x(0)-e|_{\infty}$
$\displaystyle~{}+\frac{\epsilon||G||_{\infty,*}||w_{t}||_{*}}{\eta_{1}\eta_{2}}+\sigma,$
(28) $\displaystyle|y(t)|_{*}\leq$
$\displaystyle\exp\left(-\frac{\eta_{2}t}{\epsilon}\right)|y(0)|_{*}+\frac{\epsilon||w_{t}||_{*}}{\eta_{2}}+\sigma,$
(29)

for all $t\geq 0$. The parameters $\eta_{1}$ and $\eta_{2}$ are specified in
Condition (9) and Lemma 2, respectively.

###### Proof.

For the given $(\delta,\sigma)$ and $\bar{\epsilon}>0$, let $\delta_{x}>0$ be
a constant such that

$\delta_{x}>\delta+\frac{\bar{\epsilon}||G||_{\infty,*}\delta}{\eta_{1}\eta_{2}}+\sigma.$
(30)

Note that $\delta_{x}>\delta$. By $|x(0)-e|_{\infty}\leq\delta$ and the fact
that $x(t)$ is continuous, there exists a time $T>0$ such that

$\displaystyle|x(t)-e|_{\infty}<\delta_{x},\text{ for all }t\in[0,T),$
$\displaystyle|x(T)-e|_{\infty}\geq\delta_{x}.$ (31)

Pick $\epsilon\in(0,\bar{\epsilon}]$, using Lemma 2 and Proposition 1 on
$t\in[0,T)$ yields

$|y(t)|_{*}\leq\exp\left(-\frac{\eta_{2}t}{\epsilon}\right)|y(0)|_{*}+\frac{||u_{t}||_{*}}{\eta_{2}},\text{
for all }t\in[0,T),$

where $u(t):=-\epsilon\frac{\partial h}{\partial x}f(x)+\epsilon w.$ Since
$x(t)$ is bounded on $t\in[0,T)$ and the mappings $h$ and $f$ are $C^{1}$
w.r.t. $x$, there exists a class $\mathcal{K}$ function
$\gamma:\mathbb{R}_{\geq 0}\to\mathbb{R}_{\geq 0}$ such that

$\left|-\epsilon\frac{\partial h}{\partial
x}f(x)\right|_{*}\leq\gamma(\epsilon),\text{ for all }t\in[0,T).$

Hence,

$|y(t)|_{*}\leq\exp\left(-\frac{\eta_{2}t}{\epsilon}\right)|y(0)|_{*}+\frac{\gamma(\epsilon)+\epsilon||w_{t}||_{*}}{\eta_{2}},$

for all $t\in[0,T)$. This implies that

$\|y_{t}\|_{*}\leq\delta+\frac{\epsilon\delta+\gamma(\epsilon)}{\eta_{2}},\text{
for all }t\in[0,T).$ (32)

That is, both $x(t)$ and $y(t)$ are bounded on $t\in[0,T)$.

Pick $\rho\in[0,T)$. Using Lemma 1 and Proposition 1, together with the fact
$f(e)=0$, yields

$\begin{split}|x(t)-e|_{\infty}&\leq\exp(-\eta_{1}t)|x(0)-e|_{\infty}\\\
&+||G||_{\infty,*}\int_{0}^{\rho}\exp(-\eta_{1}(t-s))|y(s)|_{*}\mathop{}\\!\mathrm{d}s\\\
&+||G||_{\infty,*}\int_{\rho}^{t}\exp(-\eta_{1}(t-s))|y(s)|_{*}\mathop{}\\!\mathrm{d}s\\\
&\leq\exp(-\eta_{1}t)|x(0)|_{\infty}+\rho||G||_{\infty,*}||y_{\rho}||_{*}\\\
&+\frac{||G||_{\infty,*}||y^{\rho}_{t}||_{*}}{\eta_{1}},\text{ for all
}t\in[0,T).\end{split}$ (33)

Note that $||y_{\rho}||_{*}\leq||y_{t}||_{*}$. Let $\rho$ be small enough such
that

$\rho||G||_{\infty,*}\left(\delta+\frac{\bar{\epsilon}\delta+\gamma(\bar{\epsilon})}{\eta_{2}}\right)\leq\frac{\sigma}{2}.$

Then, Eqs. (32) and (33) lead to

$\begin{split}&|x(t)-e|_{\infty}\leq\exp(-\eta_{1}t)|x(0)-e|_{\infty}+\frac{||G||_{\infty,*}||y_{t}^{\rho}||^{*}}{\eta_{1}}+\frac{\sigma}{2}\\\
&\leq\exp(-\eta_{1}t)|x(0)|_{\infty}+\frac{\sigma}{2}\\\
&+\frac{||G||_{\infty,*}}{\eta_{1}}\bigg{(}\exp\left(-\frac{\eta_{2}\rho}{\epsilon}\right)|y(0)|_{*}+\frac{\gamma(\epsilon)+\epsilon||w_{t}||_{*}}{\eta_{2}}\bigg{)},\end{split}$

for all $t\in[0,T)$. Note that the terms
$\exp\left(-\frac{\eta_{2}\rho}{\epsilon}\right)$ and $\gamma(\epsilon)$ goes
to zero as $\epsilon$ goes to zero. Therefore, there exists a small enough
$\epsilon_{1}>0$, such that if
$\epsilon\leq\min\\{\epsilon_{1},\bar{\epsilon}\\}$, then

$|x(t)-e|_{\infty}\leq\exp(-\eta_{1}t)|x(0)-e|_{\infty}+\frac{\epsilon||G||_{\infty,*}}{\eta_{1}\eta_{2}}||w_{t}||_{*}+\sigma,$

for all $t\in[0,T)$. From the definition of $\delta_{x}$, we have that
$|x(t)-e|_{\infty}<\delta_{x}$ for all $t\in[0,T)$. By the continuity of
$x(t)$, there exists $T^{\prime}>T$ such that
$|x(T^{\prime})-e|_{\infty}<\delta_{x}$, and this contradicts with (31).
Hence, we have $T=\infty$, and (2) holds for all $t\geq 0$. Finally, pick
$\epsilon_{2}$ such that $\gamma(\epsilon_{2})/\eta_{2}\leq\sigma$ and let
$\epsilon^{*}=\min\\{\epsilon_{1},\epsilon_{2},\bar{\epsilon}\\}$, then (2)
and (29) hold for all $t\geq 0$. $~{}~{}~{}~{}\square$

######  Remark 3.

The above proof also provides a procedure for selecting $\epsilon$
(equivalently, $k_{\operatorname{min}}$), which essentially depends on
$\delta$ and $\sigma$. Specifically, a larger $\delta$ (or a smaller $\sigma$)
requires a smaller $\epsilon$.

Similar to the formulation of $\mathcal{L}_{\infty}$ string stability (see
Def. 1 in Section 5), we seek to derive a relation between
$\max_{i}\\{||x_{i}(t)-e_{i}||\\}$ and $||w||_{\infty}$ in the next result. It
is directly follows from Theorem 2 and the definition of $|\cdot|_{*}$ in
(20), and it shows that the effect of $||w||_{\infty}$ on
$\max_{i}\\{||x_{i}(t)-e_{i}||\\}$ is scaled by
$O(\sqrt{\left\lceil\frac{n}{r}\right\rceil})$ for a fixed $n$. That is, the
effect of disturbance propagation can be reduced by increasing the
communication range $r$.

######  Corollary 1.

Consider the $(x,y)$-system (16)-(17) where the mappings $d_{i}$’s satisfy
Conditions (9)-(10), and the parameters $\bar{k}_{i}$’s satisfy Condition (21)
for some $\bar{\epsilon}>0$. Assume that $w(t)$ is essentially bounded, then
for any bounded initial conditions $x(0)$ and $y(0)$, there exists
$\epsilon^{*}>0$ such that if $\epsilon\in(0,\epsilon^{*}]$, then

$\displaystyle||x_{i}-e_{i}||\leq|x(0)-e_{i}|_{\infty}+$
$\displaystyle\epsilon
c\sqrt{\left\lceil\frac{n}{r}\right\rceil}||w||_{\infty}+\sigma,$
$\displaystyle i=1,\dots,n,$ (34)

for some $c,\sigma>0$.

###### Proof.

By Theorem 2 and (27), a small enough $\epsilon$ ensures that

$\displaystyle||x_{i}-e_{i}||\leq|x(0)-e|_{\infty}+\frac{2\epsilon||w||_{*}}{\eta_{1}\eta_{2}}+\sigma,$

for all $i=1,\dots,n$. That is, the effect of the disturbance $w(t)$ on every
$x_{i}$ is related to the norm $|\cdot|_{*}$ defined in (20). Let $q_{i}>0$,
$i=1,\dots,m$, be the $i$th diagonal entry of $D$, and define
$\bar{q}:=\max_{i}\\{q_{i}\\}$. Then for $a\in\mathbb{R}^{n}$ with the
decomposition as in Lemma 2,

$\displaystyle|a|_{*}=\left(\sum\limits_{i=1}^{m}\left(q_{i}|a^{i}|_{\infty}\right)^{2}\right)^{\frac{1}{2}}\leq\sqrt{\left\lceil\frac{n}{r}\right\rceil}\bar{q}|a|_{\infty}.$

This implies that
$||w||_{*}\leq\sqrt{\left\lceil\frac{n}{r}\right\rceil}\bar{q}||w||_{\infty}$.
Therefore, Eq. (34) holds with $c:=\frac{2\bar{q}}{\eta_{1}\eta_{2}}$.
$~{}~{}~{}~{}\square$

As shown in Example 1, the parameter $c$ in (34) also depends on $n$ and $r$
implicitly. The exact relation between $c$ and $n$ (or $r$) is not the focus
of our current work. However, for a fixed $n$, we can always find a maximal
value of $c$ such that (34) holds for all $r\in\\{1,\dots,n\\}$.

Note from (34) that a $\mathcal{L}_{\infty}$ string stability-like result is
achieved if $r=n$. Inspired by this, the next section proposes a control
design ensuring $\mathcal{L}_{\infty}$ string stability of the platooning
system (7). Specifically, it is related to an analogue of Theorem 2 but with
$\sigma=0$.

##  5 A $\mathcal{L}_{\infty}$ string stability result

This section considers a simplified version of the control protocol (12) with
$r=n$, i.e.,

$\begin{split}\frac{1}{m_{i}}u_{i}&=-k_{i}(v_{i}-v_{di})+\frac{\partial}{\partial
x_{i}}d_{i}(x_{i},x_{i+1})(v_{i-1}-v_{i})\\\ &~{}+\frac{\partial}{\partial
x_{i+1}}d_{i}(x_{i},x_{i+1})(v_{i}-v_{i+1}),~{}i=1,\dots,n.\end{split}$ (35)

Comparing (35) to (12) with $r=n$, note that the term
$\sum_{j=1}^{n-1}d_{i-j}$ (i.e., $\sum_{j=1}^{i-1}d_{j}$) is dropped. As a
result, each local controller only depend the information of tow nearest
neighbours and the velocity of the leading vehicle, i.e., $v_{0}(t)$. Compared
to the control design given by Monteil et al. (2019), the relative positions
w.r.t. the leading vehicle, i.e., $p_{0}(t)-p_{i}(t)$, $i=1,\dots,n$, are not
required.

The subsequent analysis show that the control law (35) ensures
$\mathcal{L}_{\infty}$ string stability of the platoon, which is formally
defined below.

######  Definition 1.

Consider the platooning system (7) and define the lumped error state

$\xi:=\begin{bmatrix}(x-e)^{T}&v_{1}-v_{0}&\cdots&v_{n}-v_{0}\end{bmatrix}^{T}.$
(36)

The system (7) is called $\mathcal{L}_{\infty}$ string stable if there exist
class $\mathcal{K}$ functions $\alpha_{i}$, $i=1,2$, such that for any initial
condition $\xi(0)\in\mathbb{R}^{2n}$, and any essentially bounded
$\dot{v}_{0}(t)$ and $\theta(t)$,

$\|x_{i}-e_{i}\|\leq\alpha_{1}(|\xi(0)|_{\infty})+\alpha_{2}(\|w\|_{\infty}),$
(37)

for all $i\in\\{1,\dots,n\\}$, and all $n\in\mathbb{N}_{+}$.

######  Remark 4.

If $\theta_{i}(t)\equiv 0$ for all $i$, then Definition 1 recovers the
$\mathcal{L}_{\infty}$ string stability given in (Ploeg et al., 2013, Def. 1).
Note that Condition (37) is independent of the length of the platoon.

Similar to the derivation of (16)-(17), the closed-loop system of (7) with the
control laws (35) is obtained as

$\displaystyle\dot{x}$ $\displaystyle=f_{1}(x)+Gz,$ (38)
$\displaystyle\frac{\mathop{}\\!\mathrm{d}z}{\mathop{}\\!\mathrm{d}\tau}$
$\displaystyle=-\operatorname{diag}(\bar{k}_{i})z+\epsilon w,$ (39)

where

$f_{1}(x):=\begin{bmatrix}-d_{1}(x_{1},x_{2})\\\
d_{1}(x_{1},x_{2})-d_{2}(x_{2},x_{3})\\\ \vdots\\\
d_{n-1}(x_{n-1},x_{n})-d_{n}(x_{n})\end{bmatrix}.$

######  Theorem 3.

Consider the platooning system (7) with the distributed controllers in (35),
i.e., the system (38)-(39). Pick mappings $d_{i}$’s such that Condition (10)
holds, and there exists $\eta>0$ such that

$\begin{split}&\frac{\partial d_{i}}{\partial x_{i}}>0,~{}\frac{\partial
d_{i}}{\partial x_{i+1}}<0,\\\ &\frac{\partial d_{1}}{\partial
x_{1}}+\frac{\partial d_{1}}{\partial x_{2}}\geq\eta,\\\ &\frac{\partial
d_{i}}{\partial x_{i}}+\frac{\partial d_{i}}{\partial
x_{i+1}}\leq\frac{\partial d_{i+1}}{\partial x_{i+1}}+\frac{\partial
d_{i+1}}{\partial x_{i+2}}-\eta,\end{split}$ (40)

for all $x\in\mathbb{R}^{n}$ and all $i\in\\{1,\dots,n\\}$, Then the
platooning system (7) is $\mathcal{L}_{\infty}$ string stable.

###### Proof.

Define $J_{1}(x):=\frac{\partial}{\partial x}f_{1}(x)$. Condition (40) ensures
that $\mu_{\infty}(J_{1}(x))\leq-\eta$ for all $x\in\mathbb{R}^{n}$.
Additionally, the Jacobian of the $z$-system (39) is
$-\operatorname{diag}(\bar{k}_{i})$, and
$\mu_{\infty}(-\operatorname{diag}(\bar{k}_{i}))=-1$.

Note that the $z$-system (39) is independent of $x$. That is, the
$(x,z)$-system (38)-(39) is a series connection of two contractive systems.
Recall that such a system is also contractive (see e.g., Ofir et al. (2022b)).
Specifically, this implies that $x(t)$ and $z(t)$ are uniformly bounded if
$w(t)\in\mathcal{L}_{\infty}$. Using Proposition 1 leads to

$\displaystyle|x(t)-e|_{\infty}\leq\exp(-\eta_{1}t)|x(0)|_{\infty}+\frac{||G||_{\infty}}{\eta}||z_{t}||_{\infty},$
$\displaystyle|z(t)|_{\infty}\leq\exp\left(-\frac{t}{\epsilon}\right)|z(0)|_{\infty}+\epsilon||w_{t}||_{\infty},\text{
for all }t\geq 0.$

This implies that $||z||_{\infty}\leq|z(0)|_{\infty}+\epsilon||w||_{\infty}.$
Hence,

$\begin{split}||x-e||_{\infty}&\leq|x(0)-e|_{\infty}+\frac{||G||_{\infty}}{\eta}(|z(0)|_{\infty}+\epsilon||w||_{\infty}).\end{split}$
(41)

Note that $||x_{i}-e_{i}||\leq||x-e||_{\infty}$ for all $i$. Recall that
$z_{i}:=v_{i}-d_{i}(x_{i},x_{i+1})-v_{0}$, and $d_{i}(e_{i},e_{i+1})=0$. By
Condition (10), the term
$|x(0)|_{\infty}+\frac{||G||_{\infty}|z(0)|_{\infty}}{\eta}$ can be viewed as
a weighted vector norm for $\xi(0)$. That is, there exists a constant
$c^{\prime}$ such that

$|x(0)-e|_{\infty}+\frac{||G||_{\infty}|z(0)|_{\infty}}{\eta}\leq
c^{\prime}|\xi(0)|_{\infty}.$ (42)

Combining (41) and (42) proves that Condition (37) in Definition 1 holds and
this completes the proof. $~{}~{}~{}~{}\square$

######  Remark 5.

Note that Condition (40) is less easier to be satisfied than Condition (9).
This is due to that the former is related to a heterogeneous control design,
that is, letting each vehicle have different local controller. For example,
let $d_{i}$’s be the linear form in (11). Condition (40) translates to

$\ell_{1}^{p}>\ell_{1}^{f},~{}~{}\ell_{i}^{p}-\ell_{i}^{f}<\ell_{i+1}^{p}-\ell_{i+1}^{f}.$

The above condition is satisfied if we fix $\ell_{i}^{p}$ [$\ell_{i}^{f}$] as
the same constant for all $i$, and $\ell_{i}^{f}$ [$\ell_{1}^{p}$] as a
decreasing [increasing] sequence w.r.t. $i$. Intuitively speaking, this means
that the vehicles at the tail [head] part should have a weaker [stronger]
“connection” with the follower [predecessor].

##  6 Numerical validation

The results in Section 4 are illustrated via simulations for a platooning
system consisting of $11$ vehicles using the control protocol (12) with

$d_{i}(x_{i},x_{i+1})=\ell_{i}g_{i}(x_{i},x_{i+1})+b_{i}(x_{i}-e_{i}),~{}i=1,\dots,10.\\\
$

where
$g_{i}(x_{i},x_{i+1})=\tanh(\ell_{i}^{p}(x_{i}-e_{i})-\ell_{i}^{f}(x_{i+1}-e_{i+1})).$

We consider three cases with communication range $r=1,3$, and $10$. We use the
same settings of control parameters for the three cases for a comparison
study. Specifically, the control gains $\ell_{i},\ell_{i}^{p}$,
$\ell_{i}^{f}$, and $b_{i}$ are selected such that Conditions (9)-(10) are
satisfied. Specifically, we pick
$\ell_{i}=0.5,\ell_{i}^{p}=0.18,~{}\ell_{i}^{f}=0.18,~{}b_{i}=0.1,$ for all
$i=1,\dots,10.$ The control gains $k_{i}$’s in (12) are chosen based on Lemma
2 and Theorem 2, and we set $k_{1}=\cdots=k_{10}=5.$

The simulation results with $r=1,3$, and $10$ are illustrated in Fig. 1, where
we set $e_{1}=\cdots=e_{10}=10\operatorname{m}$. The time-varying velocity of
the leading vehicle is given as

$v_{0}(t)=\begin{cases}15,&t\in[0,5),\\\ 5+2t,&t\in[5,15),\\\
35,&t\in[15,25),\\\ 85-2t,&t\in[25,35),\\\ 15,&t\in[35,45),\\\
82.5-1.5t,&t\in[45,55),\\\ 0,&t\in[55,65),\\\ -97.5+1.5t,&t\in[65,75),\\\
15,&t\in[75,\infty).\end{cases}$

The external disturbances are

$\frac{\theta_{i}}{m_{i}}=c_{i}\exp(-0.02t)\sin(t)(\text{or
}\cos(t)),~{}i=1,\dots,10,$

where $c_{i}\in[-3,3]$. Fig. 1 shows that the bounds for
$\max_{i}\\{||x_{i}-e_{i}||\\}$ is reduced by increasing the communication
range $r$. It is also observed that the effect of $\dot{v}_{0}$ attenuates
faster along the vehicular string with a larger $r$.

Figure 1: Position deviations, i.e., $x_{i}(t)-e_{i}$, $i=1,\dots,10$, of the
platooning systems with communication range $r=1,3,10$. (Blue line:
$x_{1}(t)-e_{1}$, red line: $x_{10}(t)-e_{10}$.)

##  7 Conclusions

This paper presents a contraction theory framework for analyzing disturbance
propagation in vehicular platoons with different communication ranges. The
disturbances under concern include both the external disturbances acting on
each vehicle and the disturbance caused by the acceleration of the leading
vehicle. We show that the effect of disturbances on the maximal overshoot of
spacing errors is scaled by $O(\sqrt{\left\lceil\frac{n}{r}\right\rceil})$ for
a vehicular platoon with fixed length $n+1$ and variable communication range
$r$. We also propose a heterogeneous control protocol which only rely on
communications between consecutive vehicles and broadcasting the velocity of
the leading vehicle such that $\mathcal{L}_{\infty}$ string stability is
achieved. We believe the proposed framework can also be utilized to analyze
other generic control protocols, for example, with mixed communication ranges.

## References

  * Aminzare and Sontag [2014] Z. Aminzare and E. D. Sontag.  Contraction methods for nonlinear systems: A brief introduction and some open problems.  In _Proc. 53rd IEEE Conf. on Decision and Control_ , pages 3835–3847, Los Angeles, CA, 2014.
  * Barooah and Hespanha [2005] P. Barooah and J. P. Hespanha.  Error amplification and disturbance propagation in vehicle strings with decentralized linear control.  In _Proc. 44th IEEE Conf. on Decision and Control_ , pages 4964–4969. IEEE, 2005.
  * Besselink and Johansson [2017] B. Besselink and K. H. Johansson.  String stability and a delay-based spacing policy for vehicle platoons subject to disturbances.  _IEEE Trans. Automat. Control_ , 62(9):4376–4391, 2017.
  * Christofides and Teel [1996] P. D. Christofides and A. R. Teel.  Singular perturbations and input-to-state stability.  _IEEE Trans. Automat. Control_ , 41(11):1645–1650, 1996.
  * Desoer and Haneda [1972] C. Desoer and H. Haneda.  The measure of a matrix as a tool to analyze computer algorithms for circuit analysis.  _IEEE Trans. Circuit Theory_ , 19(5):480–486, 1972.  10.1109/TCT.1972.1083507.
  * Farnam and Sarlette [2019] A. Farnam and A. Sarlette.  Toward a comprehensive impossibility result for string stability.  _IEEE Trans. Automat. Control_ , 65(4):1652–1659, 2019.
  * Feng et al. [2019] S. Feng, Y. Zhang, S. E. Li, Z. Cao, H. X. Liu, and L. Li.  String stability for vehicular platoon control: Definitions and analysis methods.  _Annu Rev Control_ , 47:81–97, 2019.
  * Herman et al. [2016] I. Herman, D. Martinec, Z. Hurák, and M. Sebek.  Scaling in bidirectional platoons with dynamic controllers and proportional asymmetry.  _IEEE Trans. Automat. Control_ , 62(4):2034–2040, 2016.
  * Herman et al. [2017] I. Herman, S. Knorn, and A. Ahlén.  Disturbance scaling in bidirectional vehicle platoons with different asymmetry in position and velocity coupling.  _Automatica_ , 82:13–20, 2017.
  * Khalil [2002] H. K. Khalil.  _Nonlinear Systems_.  Prentice hall, Upper Saddle River, NJ, third edition, 2002.
  * Li et al. [2019] S. E. Li, X. Qin, Y. Zheng, J. Wang, K. Li, and H. Zhang.  Distributed platoon control under topologies with complex eigenvalues: Stability analysis and controller synthesis.  _IEEE Trans. Control Syst. Technol_ , 27(1):206–220, 2019.
  * Mesbahi and Egerstedt [2010] M. Mesbahi and M. Egerstedt.  Graph theoretic methods in multiagent networks.  In _Graph Theoretic Methods in Multiagent Networks_. Princeton University Press, 2010.
  * Monteil et al. [2019] J. Monteil, G. Russo, and R. Shorten.  On $\mathcal{L}_{\infty}$ string stability of nonlinear bidirectional asymmetric heterogeneous platoon systems.  _Automatica_ , 105:198–205, 2019.
  * Ofir et al. [2022a] R. Ofir, F. Bullo, and M. Margaliot.  Minimum effort decentralized control design for contracting network systems.  _IEEE Control Systems Letters_ , 6:2731–2736, 2022a.
  * Ofir et al. [2022b] R. Ofir, M. Margaliot, Y. Levron, and J.-J. Slotine.  A sufficient condition for k-contraction of the series connection of two systems.  _IEEE Trans. Automat. Control_ , 67(9):4994–5001, 2022b.
  * Peters et al. [2014] A. A. Peters, R. H. Middleton, and O. Mason.  Leader tracking in homogeneous vehicle platoons with broadcast delays.  _Automatica_ , 50(1):64–74, 2014.
  * Peters et al. [2016] A. A. Peters, R. H. Middleton, and O. Mason.  Cyclic interconnection for formation control of 1-d vehicle strings.  _European J. Control_ , 27:36–44, 2016.
  * Ploeg et al. [2013] J. Ploeg, N. Van De Wouw, and H. Nijmeijer.  $\mathcal{L}_{p}$ string stability of cascaded systems: Application to vehicle platooning.  _IEEE Trans. Control Syst. Technol_ , 22(2):786–793, 2013.
  * Russo et al. [2010] G. Russo, M. di Bernardo, and E. D. Sontag.  Global entrainment of transcriptional systems to periodic inputs.  _PLOS Computational Biology_ , 6:e1000739, 2010.
  * Russo et al. [2013] G. Russo, M. di Bernardo, and E. D. Sontag.  A contraction approach to the hierarchical analysis and design of networked systems.  _IEEE Trans. Automat. Control_ , 58(5):1328–1331, 2013.
  * Seiler et al. [2004] P. Seiler, A. Pant, and K. Hedrick.  Disturbance propagation in vehicle strings.  _IEEE Trans. Automat. Control_ , 49(10):1835–1842, 2004.
  * Silva et al. [2021] G. F. Silva, A. Donaire, A. McFadyen, and J. J. Ford.  String stable integral control design for vehicle platoons with disturbances.  _Automatica_ , 127:109542, 2021.
  * Sontag et al. [1989] E. D. Sontag et al.  Smooth stabilization implies coprime factorization.  _IEEE Trans. Automat. Control_ , 34(4):435–443, 1989.
  * Stankovic et al. [2000] S. S. Stankovic, M. J. Stanojevic, and D. D. Siljak.  Decentralized overlapping control of a platoon of vehicles.  _IEEE Trans. Control Syst. Technol_ , 8(5):816–832, 2000.
  * Ström [1975] T. Ström.  On logarithmic norms.  _SIAM J. Numerical Analysis_ , 12(5):741–753, 1975.
  * Stüdli et al. [2018] S. Stüdli, M. M. Seron, and R. H. Middleton.  Vehicular platoons in cyclic interconnections.  _Automatica_ , 94:283–293, 2018.
  * Swaroop and Hedrick [1996] D. Swaroop and J. K. Hedrick.  String stability of interconnected systems.  _IEEE Trans. Automat. Control_ , 41(3):349–357, 1996.
  * Varaiya [1993] P. Varaiya.  Smart cars on smart roads: problems of control.  _IEEE Trans. Automat. Control_ , 38(2):195–207, 1993.
  * Yadlapalli et al. [2006] S. K. Yadlapalli, S. Darbha, and K. R. Rajagopal.  Information flow and its relation to stability of the motion of vehicles in a rigid formation.  _IEEE Trans. Automat. Control_ , 51(8):1315–1319, 2006.
  * Yueh [2005] W.-C. Yueh.  Eigenvalues of several tridiagonal matrices.  _Applied Mathematics E-Notes_ , 5:66–74, 2005.
  * Zheng et al. [2015] Y. Zheng, S. E. Li, K. Li, and L.-Y. Wang.  Stability margin improvement of vehicular platoon considering undirected topology and asymmetric control.  _IEEE Trans. Control Syst. Technol_ , 24(4):1253–1265, 2015.