\documentclass[12pt]{article} \usepackage{e-jc} \specs{P15}{19(3)}{2012} \usepackage{eepic,amsmath,amsmath} \newtheorem{theorem}{Theorem} \newtheorem{lemma}{Lemma} \newtheorem{corollary}{Corollary} \newtheorem{definition}{Definition} \newtheorem{proposition}{Proposition} \newtheorem{result}{Result} \newtheorem{conjecture}{Conjecture} \newcommand{\con}{\rightarrow} \newcommand{\cd}{\mbox{$\stackrel{d}{\rightarrow}$}} \newcommand{\ep}{\varepsilon} \newcommand\be{\begin{equation}} \newcommand\ee{\end{equation}} \newcommand\ba{\begin{array}{lcr}} \newcommand\ea{\end{array}} \newcommand\bea{\begin{eqnarray*}} \newcommand\eea{\end{eqnarray*}} \newcommand\QED{\ifhmode\allowbreak\else\nobreak\fi\quad\nobreak$\Box$\medbreak} \newcommand{\proofstart}{\par\noindent\bf Proof:\rm\enspace} \newcommand{\proofend}{\hfill\eod\med} \newenvironment{proof}{\proofstart}{\proofend} \def\noi{\noindent} \def\parsec{\par\noindent} \def\big{\bigskip\parsec} \def\med{\medskip\parsec} %\def\basswidth{3.0in} %\def\gutter{0.1in} \def\noi{\noindent} \def\eod{\vrule height 6pt width 5pt depth 0pt} \def\parsec{\par\noindent} \def\big{\bigskip\parsec} \def\med{\medskip\parsec} \def\pr{{\rm Pr}} \def\l{\ell} \def\A{{\cal A}} \def\kk{\kappa} \def\W{{\cal W}} \def\X{\{X_k\}} \def\Wk{{\A_{k}}} \def\eps{\varepsilon} \def\B{{\cal B}} \def\TB{\widetilde{\B}} \def\F{{\cal F}} \def\S{{\cal S}} \def\TG{\widetilde{G}} \def\P{{\cal P}} \def\HC{{\cal H}} \def\tx{\hat{x}} \def\HX{\hat{X}} \def\tb{\widetilde{B}} \def\tw{\widetilde{\W}} \def\TX{\widetilde{X}} \def\pmin{P_{\min}} \def\lf{\lfloor} \def\rf{\rfloor} \def\lb{L_n^{(b)}} \def\hb{H_n^{(b)}} \def\sb{s_n^{(b)}} \def\O{{\cal O}} \def\var{{\rm Var~}} \def\cov{{\rm Cov}} \def\J{{\cal J}} \def\ow{\overline{w}} \def\ta{{\widetilde{a}}} \def\real{{I \!\! R}} \def\h{\hat} \def\D{{\cal D}} \def\Var{{\rm Var}~} \def\la{\langle} \def\ra{\rangle} \def\Var{{\bf Var~}} \def\E{{\bf E}} \def\vpi{\mbox{\boldmath $\pi$}} \def\wt{\widetilde} \def\t{\theta} \def\dd{\delta} \def\TX{\wt{X}} \def\TV{\wt{V}} \def\C{{\cal C}} \def\d{\delta} \def\bb{\beta} \def\OG{O_n^{{\rm gr}}} \def\OP{O_n^{{\rm opt}}} \def\G{{\cal G}} \def\r{\bar{r}} \def\M{{\cal M}} \def\T{{\cal T}} \def\U{{\cal U}} \def\L{{\cal L}} \def\R{{\cal R}} \def\om{\omega} \def\tb{\tilde{b}} \def\F{{\cal F}} \def\Q{{\bf Q}} \def\V{{\bf V}} \def\N{{\cal N}} \def\K{{\cal K}} \def\A{{\cal A}} \def\M{{\cal M}} \def\Z{{\cal Z}} \def\k{{\bf k}} \def\kk{{\bf k'}} \def\skk{{\k,\kk}} \def\skn{{\k,[0]}} \def\skm{{\k-\kk}} \def\det{\hbox{\rm det}} \def\TG{\tilde{G}} \def\j{{\bf i}} \def\ss{\sigma} \def\eps{\epsilon} \def\z{{\bf z}} \def\x{{\bf x}} \def\v{{\bf v}} \def\kb{{\bf k}} \def\y{{\bf y}} \def\t{{\bf t}} \title{\bf On a Recurrence Arising\\ in Graph Compression\thanks{This work was supported in part by the NSF Science and Technology Center for Science of Information Grant CCF-0939370, NSF Grant CCF-0830140, AFOSR Grant FA8655-11-1-3076, NSA Grants H98230-11-1-0141 and H98230-11-1-0184, and the MNSW grant N206 369739.}} % input author, affilliation, address and support information as follows; % the address should include the country, and does not have to include % the street address \author{Yongwook Choi\\ \small J.~Craig Venter Institute\\[-0.8ex] \small Rockville, MD, USA 20850\\ \small\tt ychoi@jcvi.org \and Charles Knessl\\ \small Dept. Math. Stat. \& Comp. Sci.\\[-0.8ex] \small University of Illinois at Chicago\\[-0.8ex] \small Chicago, IL, USA 60607\\ \small\tt knessl@uic.edu \and Wojciech Szpankowski\thanks{Also, Visiting Professor at ETI, Gdansk University of Technology, Poland.}\\ \small Dept. Computer Science\\[-0.8ex] \small Purdue University\\[-0.8ex] \small West Lafayette, IN, USA 47907\\ \small\tt spa@cs.purdue.edu } \date{\dateline{November 16, 2011}{July 11, 2012}{Aug 9, 2012}\\ \small Mathematics Subject Classification: 05C05, 68W40} \begin{document} \maketitle \begin{abstract} In a recently proposed graphical compression algorithm by Choi and Szpankowski (2012), the following tree arose in the course of the analysis. The root contains $n$ balls that are consequently distributed between two subtrees according to a simple rule: In each step, all balls independently move down to the left subtree (say with probability $p$) or the right subtree (with probability $1-p$). A new node is created as long as there is at least one ball in that node. Furthermore, a nonnegative integer $d$ is given, and at level $d$ or greater one ball is removed from the leftmost node before the balls move down to the next level. These steps are repeated until all balls are removed (i.e., after $n+d$ steps). Observe that when $d=\infty$ the above tree can be modeled as a {\it trie} that stores $n$ independent sequences generated by a binary memoryless source with parameter $p$. Therefore, we coin the name $(n,d)$-tries for the tree just described, and to which we often refer simply as $d$-tries. We study here in detail the path length, and show how much the path length of such a $d$-trie differs from that of regular tries. We use methods of analytic algorithmics, from Mellin transforms to analytic poissonization. \end{abstract} \section{Introduction} In \cite{isit09-graph} an algorithm was described to compress the {\it structure} of a graph. The main idea behind the algorithm is quite simple: First, a vertex of a graph, say $v_1$, is selected and the number of neighbors of $v_1$ is stored in a binary string. Then the remaining $n-1$ vertices are partitioned into two sets: the neighbors of $v_1$ and the non-neighbors of $v_1$. This process continues by selecting randomly a vertex, say $v_2$, from the neighbors of $v_1$ and storing two numbers: the number of neighbors of $v_2$ among each of the above two sets. Then the remaining $n-2$ vertices are partitioned into four sets: the neighbors of both $v_1$ and $v_2$, the neighbors of $v_1$ that are non-neighbors of $v_2$, the non-neighbors of $v_1$ that are neighbors of $v_2$, and the non-neighbors of both $v_1$ and $v_2$. This procedure continues until all vertices are processed. In the Erd\H{o}s-R\'enyi model, a random graph has any pair of vertices connected by an edge with probability $p$. It is proved in \cite{isit09-graph} that for large $n$ our algorithm optimally compresses any graph generated by the Erd\H{o}s-R\'enyi model (and, in fact, it works well in practice even for graphs not generated by the Erd\H{o}s-R\'enyi model). To establish this asymptotic optimality result, an interesting tree was used in the construction, that we describe next. The root of such a tree contains $n$ balls (vertices of the underlying graph) that are consequently distributed between two subtrees according to a simple rule: In each step, all balls independently move down to the left subtree (say with probability $p$) or the right subtree (with probability $q=1-p$), and a new node is created as long as there is at least one ball in that node. Finally, a non-negative integer $d$ is given so that at level $d$ or greater one ball is removed from the leftmost node before the balls move down to the next level. These steps are repeated until all balls are removed (i.e., after $n$+$d$ steps). Of interest are such tree parameters as the depth, path length (sum of all depths), size, and so forth. This is illustrated in Figure~\ref{fig:trie} in which the deleted ball is shown next to the node from where it was removed. \begin{figure}[tbp] \setlength{\unitlength}{1.5mm} \center \begin{picture}(55,50)(10,20) %\multiput(0,0)(0,7){10}{\dashbox{1.5}(80,5){}}% \begin{scriptsize} \put(50,51){\circle{6}}% 3 \put(49,51){\circle*{1}}% \put(51,51){\circle*{1}}% \put(27.5,51){\circle{6}}% 3 \put(27.5,51){\circle*{1}}% \put(31,44){\circle{6}}% 4 \put(31,44){\circle*{1}}% \put(22.5,44){\circle{6}}% 4 \put(18,44){\circle*{1}}% \put(34,37){\circle{6}}% 5 \put(29.5,37){\circle*{1}}% \put(18.5,51){\circle{6}}% 3 \put(18.5,51){\circle*{1}}% \put(14,51){\circle*{1}}% \put(42,58){\circle{6}}% 2 \put(41,58){\circle*{1}}% \put(43,58){\circle*{1}}% \put(23,58){\circle{6}}% 2 \put(22,57){\circle*{1}}% \put(24,57){\circle*{1}}% \put(18.5,58){\circle*{1}}% \put(23,59){\circle*{1}}% \put(32,65){\circle{6}}% 1 \put(30.5,66){\circle*{1}}% \put(32,66){\circle*{1}}% \put(33.5,66){\circle*{1}}% \put(30.5,64){\circle*{1}}% \put(32,64){\circle*{1}}% \put(33.5,64){\circle*{1}}% \put(55,44){\circle{6}}% 4 \put(55,44){\circle*{1}}% \put(45,44){\circle{6}}% 4 \put(45,44){\circle*{1}}% \put(42,37){\circle{6}}% 5 \put(42,37){\circle*{1}}% \put(45,30){\circle{6}}% 6 \put(40.5,30){\circle*{1}}% \put(60,37){\circle{6}}% 5 \put(60,37){\circle*{1}}% \put(55,30){\circle{6}}% 6 \put(55,30){\circle*{1}}% \put(59,23){\circle{6}}% 7 \put(54.5,23){\circle*{1}}% \end{scriptsize} % between circles \put(24.1,60.7){\line(10,7){5}}% 1-2 \put(34.9,64){\line(12,-7){5.8}}% 1-2 \put(25.2,56){\line(8,-7){2.4}}% 2-3 \put(20.7,56){\line(-8,-7){2.4}}% 2-3 \put(45,57){\line(9,-7){4}}% 2-3 \put(48,48.7){\line(-6,-7){1.7}}% 3-4 \put(52,48.6){\line(6,-7){1.7}}% 3-4 \put(29.2,48.6){\line(6,-7){1.5}}% 3-4 \put(20.3,48.6){\line(6,-7){1.5}}% 3-4 \put(57.2,41.8){\line(5,-7){1.5}}% 4-5 \put(42.9,41.9){\line(-4,-7){1.1}}% 4-5 \put(32.5,41.4){\line(5,-7){1}}% 4-5 \put(43.5,34.4){\line(4,-7){0.9}}% 5-6 \put(57.6,35.2){\line(-4,-7){1.4}}% 5-6 \put(56.6,27.5){\line(4,-7){1.1}}% 6-7 \end{picture} \caption{\small A $(6,1)$-trie with six balls and $d=1$, in which the deleted ball is shown next to the node where it was removed.} \label{fig:trie} \end{figure} The tree just described falls between two digital trees, namely {\it tries} and {\it digital search trees}. In fact, when $d=\infty$ the tree can be modeled as a {\it trie} that stores $n$ independent sequences generated by a binary memoryless source with parameter $p$. Hence, we coin the term $(n,d)$-trie (or simply $d$-trie) for the tree just described. In \cite{isit09-graph} lower and upper bounds were proved for parameters of interest, by using known results for tries and digital search trees \cite{drmota09,spa-book}. In this paper, we establish precise asymptotic results. In particular, we show by how much the path length of a $d$-trie differs from the path length of the corresponding regular trie. Many parameters of a $d$-trie can be described by the following two dimensional recurrence \begin{equation} \label{spa1} a(n,d) = f(n) +\sum_{k=0}^n {n\choose k} p^kq^{n-k} [a(k,d-1)+a(n-k,k+d-1)],\;\; d\geq 1, \end{equation} with $q=1-p$, and the boundary equation \begin{equation} \label{spa2} a(n+1,0) = f(n) +\sum_{k=0}^n {n\choose k} p^kq^{n-k} [a(k,0)+a(n-k,k)], \end{equation} for a known additive term $f(n)$. For example, when $f(n)=n$, then $a(n,d)$ represents the path length. Recurrence (\ref{spa2}) is equivalent to the following boundary condition $$ a(n,1) = a(n+1,0). $$ For $d=\infty$ recurrence (\ref{spa1}) becomes a traditional recurrence arising in the analysis of tries \cite{spa-book} whose solutions (exact and asymptotic) are well known. Thus, it is natural to study the difference $\tilde{a}(n,d) := a(n,\infty)-a(n,d)$, and that is our objective. In passing, we should point out that recurrence (\ref{spa2}) resembles the one used to analyze another digital search tree, known as a {\it digital search tree}. In this paper we prove, however, that a $(n,d)$-trie more closely resembles a trie, rather than a digital search tree. Our main interest lies in solving recurrence (\ref{spa1}) for $d=O(1)$. In fact, for graph compression we only need $d=0$, and we focus on this case. We shall show that the term in (\ref{spa1}) involving the sum over $a(n-k,k+d-1)$ becomes exponentially small for $n$ large and $d$ fixed. Then we shall approximate the recurrence for the excess quantity $\tilde{a}(n,d)$ by $$ \tilde{a}(n,d) = \sum_{k=0}^n {n\choose k} p^kq^{n-k} \tilde{a}(k,d-1) $$ with an appropriate initial condition. The above we can solve asymptotically using Mellin transform technique and depoissonization. In particular, for $f(n)=n$ (that is, for the path length in a $d$-trie) we prove that the excess quantity $\tilde{a}(n,d)$ becomes asymptotically, as $n\to \infty$ and $d=O(1)$, $$ \frac{1}{2h\log (1/p)} \log^2 n - \frac{d}{h}\log n + \left[\frac{1}{2h} - \frac{1}{h\log p} \left(\gamma+1+\frac{h_2}{2h} + \Psi(\log_p n) \right)\right]\log n $$ where $\Psi(\cdot)$ is the periodic function when $\log p/\log (1-p)$ is rational, and $h$ is the entropy rate. Digital trees such as tries and digital search trees have been intensively studied for the last thirty years \cite{devroye92,drmota09,fgd,fs-book,ks00,KnuthV1,louchard87,pittel85,pittel86,spa91,spa-book}. However, our two-dimensional recurrence seems to be new and harder to analyze. It somewhat resembles the profile recurrences for digital trees, which were studied for tries in \cite{park08} and digital search trees in \cite{ds11}, and which are known to also be challenging. The paper is organized as follows. In the Section~\ref{sec:Problem Statement} we precisely formulate our problem and analyze it for $f(n)=n$. Proofs are presented in Section~\ref{sec:analysis}, where we also discuss some asymptotics for $d\to \infty$. %Some numerical studies are carried out in Section~\ref{sec:Numerical Data}. \section{Problem Statement} \label{sec:Problem Statement} In this section, we first formulate some recurrences describing $(n,d)$-tries, then summarize our main results, discuss some extensions, and present numerical results. \subsection{Main Results} Let us consider a $(n,d)$-trie with $n$ balls and parameter $d\geq 0$. First, we analyze the average path length $b(n,d)$. It satisfies the following recurrence equations \begin{equation}\label{b0} b(n+1,0) = n + \sum_{k=0}^{n} {n \choose k} p^k q^{n-k} \left[b(k,0)+b(n-k,k)\right], \ \ \mbox{for $n\ge 2$}, \end{equation} and \begin{equation}\label{bd} b(n,d) = n + \sum_{k=0}^{n} {n \choose k} p^k q^{n-k} \left[b(k,d-1)+b(n-k,k+d-1)\right], \ \ \mbox{for $n\ge 2$, $d\ge 1$}. \end{equation} Recurrence (\ref{b0}) follows from the fact that starting with $n+1$ balls in the root node, and removing one ball, we are left with $n$ balls passing through the root node. The root contributes $n$ since each time a ball moves down it adds $1$ to the path length. Those $n$ balls move down to the left or the right subtrees. Let us assume $k$ balls move down to the left subtree (the other $n-k$ balls must move down to the right subtree); this occurs with probability ${n \choose k} p^k q^{n-k}$. At level one, one ball is removed from those $k$ balls in the root of the left subtree. This contributes $b(k,0)$. There will be no removal from $n-k$ balls in the right subtree until all $k$ balls in the left subtree are removed. This contributes $b(n-k,k)$. Similarly, for $d>0$ we arrive at recurrence (\ref{bd}). Here $0 0. \end{equation} Thus $\bar{C}(\cdot)$ is a periodic function of $\log_p n$ of period 1, which we can write as the Fourier series \begin{equation}\label{Cn-Fs} \bar{C}(n) = c^{(0)}(p) + \sum_{\l=-\infty,\l\neq 0}^\infty c^{(\l)}(p)e^{2\pi i \l \log_p n}. \end{equation} It again appears difficult to identify explicitly the Fourier coefficients $c^{(\l)}(p)$, but we can do this in the symmetric case $p=q=1/2$. Then we set $\sum_{n=0}^\infty C(n)z^n/n! = \F_2(z)$ as in (\ref{FJz}) and from (\ref{BnJ2}) obtain \begin{equation}\label{Cn2} C(n) = \frac{2n!}{2\pi i} \oint \frac{e^z}{z^{n-1}}\prod_{\l=1}^\infty\left(\frac{1+e^{-z2^{-\l-1}}}{2}\right)dz. \end{equation} To obtain the large $n$ behavior of the integral in (\ref{Cn2}) we first expand the integral for $z\to\infty$ and apply a depoissonization argument. Setting $\l=\log_2 z +J$ we have $2^\l=2^Jz$ and \begin{eqnarray*} &&\prod_{\l=1}^\infty\left(\frac{1+e^{-z2^{-\l-1}}}{2}\right)\\ &=& \exp\left[\sum_{\l=1}^\infty\log \left(\frac{1+e^{-z2^{-\l-1}}}{2}\right)\right]\\ &=&\exp\left[\sum_{J=1-\log_2 z}^\infty\log \left(\frac{1+e^{-2^{-J-1}}}{2}\right)\right]\\ &=&\exp\left[\sum_{J=0}^\infty\log \left(\frac{1+e^{-2^{-J-1}}}{2}\right) + \sum_{J=1}^{\log_2z -1}\log \left(\frac{1+e^{-2^{J-1}}}{2}\right)\right]\\ &\sim&\exp\left[\log(\frac{1}{2})(\log_2z-1) + \sum_{J=0}^\infty\log \left(\frac{1+e^{-2^{-J-1}}}{2}\right) + \sum_{J=1}^{\infty}\left(1+e^{-2^{J-1}}\right)\right]\\ &=&\frac{2}{z}K_* \end{eqnarray*} where $$ K_* =\prod_{J=0}^\infty \left(\frac{1+e^{-2^{-J-1}}}{2}\right) \prod_{J=1}^{\infty}\exp\left(1+e^{-2^{J-1}}\right) = 1. $$ Thus $C(n) \sim 4n!/(n-1)! = 4n$ as $n\to\infty$. This shows that $c^{(0)}(1/2) = 4$ and a more careful calculation can be used to identify the other Fourier coefficients $c^{(\l)}(1/2)$ in (\ref{Cn-Fs}) (then we would set $\l=\lfloor\log_2 z\rfloor+J = \log_2 z +J -\{\log_2 z\}$ so that $2^\l=2^Jz2^{-\{\log_2 z\}}$. We omit the details. In Table~\ref{table:cn} we consider various values of $p$ and estimate $\bar{C}(n) \approx c^{(0)}(p)$ numerically, by computing $C(n)n^{-\nu}$ from (\ref{Cn}), for large $n$. This shows that as a function of $p$, $|c^{(0)}(p)|$ is minimal when $p$ is between 0.6 and 0.7, and becomes large as either $p\to 0$ or $p\to 1$. For $p\to 0$ the oscillatory terms in (\ref{Cn-Fs}) become more numerically significant. Table~\ref{table:cn} indicates this when $p=0.1$, by giving a range of values of $C(n)n^{-\nu}$. \begin{table}[htb] \begin{center} \caption{Values of the zeroth Fourier coefficient.\label{table:cn}} \par \medskip \begin{tabular}{|c|c|} \hline % after \\: \hline or \cline{col1-col2} \cline{col3-col4} ... $p$ & $C(n)n^{-\nu}|_{n\to\infty}\approx c^{(0)}(p) $ \\ \hline 0.5 & 4 \\ 0.4 & 5.664 \\ 0.3 & 9.728\\ 0.25 & 14.03\\ 0.2 & 22.5\\ 0.1 & 98 to 105\\ & \\ 0.6 & 3.331\\ 0.7 & 3.276\\ 0.75 & 3.479\\ 0.8 & 3.903\\ 0.9 & 6.423\\ \hline \end{tabular} \end{center} \end{table} \par \medskip To justify the approximation in (\ref{Cn-sim}) we first inductively show that for all $n$ \begin{equation}\label{Cn-le} C(n) \le An^{\nu+\eps} \end{equation} for all $\eps>0$ and $A>0$. By isolating the terms in the sums in (\ref{Cn}) with $k=n$ and $k=0$ we obtain, for $n>2$, \begin{equation}\label{Cn2a} C(n) = \frac{1}{p^2+q^2-p^n-q^n}\left[\sum_{k=2}^{n-1} {n\choose k}p^kq^{n-k}C(k) + \sum_{k=1}^{n-2} {n\choose k}p^kq^{n-k} (p^2+q^2)^kC(n-k)\right]. \end{equation} Assuming inductively that (\ref{Cn-le}) holds for $C(k)$ for $k=1,2,\cdots,n-1$ we then have \begin{eqnarray*} \sum_{k=2}^{n-1} {n\choose k}p^kq^{n-k}C(k) &\le & \sum_{k=2}^{n-1} {n\choose k}p^kq^{n-k}Ak^{\nu+\eps} \\ &\le& A(np)^{\nu+\eps}. \end{eqnarray*} Using a similar estimate for the second sum in (\ref{Cn2a}) we are led to \begin{eqnarray}\label{Cn3} C(n) &\le & \frac{A}{p^2+q^2-p^n-q^n}\left[(np)^{\nu+\eps} + n^{\nu+\eps}(p(p^2+q^2)+q)^n\right] \nonumber \\ &=&An^{\nu+\eps}\left[\frac{p^2+q^2}{p^2+q^2-p^n-q^n}p^{\eps} +\frac{(p(p^2+q^2)+q)^n}{p^2+q^2-p^n-q^n} \right], \end{eqnarray} as $C(n-k)\le A(n-k)^{\nu+\eps} \le An^{\nu+\eps}$ and $p^\nu = p^2+q^2$. Since $p(p^2+q^2)+q < p+q = 1$, the second term in (\ref{Cn3}) is asymptotically negligible for $n$ large and (\ref{Cn-le}) follows by induction. We have thus obtained some exact expressions for $b(n,d)$ for small values of $n$, a general asymptotic result for $d\to\infty$ with $n=O(1)$, and then examined how this result behaves when $n$ also becomes large. However, this cannot be used to infer the behavior of $b(n,d)$ for $n\to\infty$ with $d=O(1)$, which we examine next. \subsection{Main Asymptotic Result for $b(n,d)$} We first give an intuitive derivation of the asymptotics of $b(n,d)$ for fixed $d\ge 0$ and $n\to\infty$, and in particular of $b(n,0)$. Starting from (\ref{btd}) we again argue that the second sum is negligible for $n\to\infty$ and that the first is asymptotic to $\tb(np,d-1)$ so that (\ref{btd}) becomes \begin{equation}\label{btnd2} \tb(n,d) \sim \tb(np,d-1), \ n\to\infty \end{equation} and, in particular, \begin{equation}\label{btn1} \tb(n,1) \sim \tb(np,0), \ n\to\infty \end{equation} which when added to (\ref{bt-boundary}) leads to \begin{equation}\label{btn0} \tb(n+1,0)-\tb(np,0) \sim b(n+1,\infty) - b(n,\infty). \end{equation} The right side of (\ref{btn0}) may be estimated from (\ref{bninf4}) or by (\ref{bninf3}). Using (\ref{bninf3}) we can show that term by term differentiating of the asymptotic series in (\ref{bninf4}) is permissible, and thus (\ref{btn0}) becomes, for $n\to\infty$, \begin{equation}\label{btn02} \tb(n+1,0)-\tb(np,0) = \frac{1}{h}\log n + \frac{1}{h}\left(\gamma+1+\frac{h_2}{2h}\right) + \frac{1}{h}\psi(\log_p n)+o(1), \end{equation} where $\psi(\cdot)$ is the periodic function \begin{equation}\label{psi} \psi(x) = \sum_{k=-\infty, k\neq 0}^\infty \left[1+\frac{2k\pi ir}{\log p}\right]\Gamma\left(-\frac{2k\pi ir}{\log p}\right)e^{2k\pi irx}, \end{equation} where we note that $\psi$ and $\Phi$ are related by $\psi(x)=\Phi(x)+(\log p)^{-1}\Phi'(x)$. Now (\ref{btn02}) suggests that $\tb(n,0)$ admits an asymptotic expansion of the form \begin{equation}\label{btn03} \tb(n,0) = A\log^2 n +B\log n+C+o(1), \ n\to\infty \end{equation} and then \begin{equation}\label{btn04} \tb(n+1,0)-\tb(np,0) =-2A(\log p)\log n -A\log^2p-B\log p +o(1). \end{equation} Comparing (\ref{btn02}) to (\ref{btn04}) we conclude that $A=-(2h\log p)^{-1}$ and then \begin{equation}\label{B} B= \frac{1}{2h} -\frac{1}{h\log p}\left[\gamma+1+\frac{h_2}{2h}+\psi(\log_p n)\right]. \end{equation} We have thus formally derived the result in Theorem~\ref{th:main} for $\tb(n,0)$. For any fixed $d>0$ we can extend this argument by asymptotically solving (\ref{btnd2}) by an expansion of the form \begin{equation}\label{btnd3} \tb(n,d) = A(d)\log^2 n +B(d)\log n+C(d)+o(1) \end{equation} to find from (\ref{btnd2}) that $A(d)=A(d-1)$ and $B(d)=B(d-1)+2\log pA(d-1)$. Then using (\ref{btnd3}) in (\ref{btn0}) or (\ref{btn02}) we find that $A(d)=A(0)=-(2h\log p)^{-1}$ and $B(d)-B(d-1) =2\log pA(d-1)=-h^{-1}$ so that $B(d)=B(0)-h^{-1}d$, where $B(0)=B$ is as in (\ref{B}). We proceed to provide a rigorous derivation of the theorem. Using arguments completely analogous to (\ref{Cn-le})--(\ref{Cn3}), we can inductively establish the bound \begin{equation}\label{btnd4} \tb(n,d) \le A_0 n^{\nu+\eps}(p^2+q^2)^d;\ n\ge 2, d\ge 0 \end{equation} where again $\nu$ is given by (\ref{nu}). When $n=2$ we have (exactly) $$\tb(2,d)=\left(\frac{1}{pq}-2\right)(p^2+q^2)^{d-1},$$ so (\ref{btnd4}) clearly holds. Assuming that (\ref{btnd4}) holds for all $(N,D)$ with $N+D