Computing the $p$-Spectral Radii of Uniform Hypergraphs with Applications

TL;DR

This paper introduces a novel conjugate gradient algorithm, CSRH, for efficiently computing the p-spectral radius of uniform hypergraphs, with applications in large-scale data analysis and hypergraph optimization.

Contribution

The paper develops a globally convergent algorithm, CSRH, for calculating the p-spectral radius of hypergraphs, outperforming existing methods especially on large-scale problems.

Findings

CSRH efficiently computes p-spectral radii for hypergraphs with millions of vertices.

CSRH reliably finds the global maximizer due to the semialgebraic structure.

The method effectively ranks real-world datasets based on hypergraph spectral properties.

Abstract

The $p$-spectral radius of a uniform hypergraph covers many important concepts, such as Lagrangian and spectral radius of the hypergraph, and is crucial for solving spectral extremal problems of hypergraphs. In this paper, we establish a spherically constrained maximization model and propose a first-order conjugate gradient algorithm to compute the $p$-spectral radius of a uniform hypergraph (CSRH). By the semialgebraic nature of the adjacency tensor of a uniform hypergraph, CSRH is globally convergent and obtains the global maximizer with a high probability. When computing the spectral radius of the adjacency tensor of a uniform hypergraph, CSRH stands out among existing approaches. Furthermore, CSRH is competent to calculate the $p$-spectral radius of a hypergraph with millions of vertices and to approximate the Lagrangian of a hypergraph. Finally, we show that the CSRH method is capable of ranking real-world data set based on solutions generated by the $p$-spectral radius model.

Figures

Tables

Table 1. Table 1: Z-Eigenvalues of adjacency tensors of several small hypergraphs.

Hypergraph	CSRH				SS-HOPM
	Iter.	Time(s)	Accu.	Err.	Iter.	Time(s)	Accu.	Err.
$G_1$	13593	3.35	1.00	$5.44\times {10}^{-16}$	2668	4.89	1.00	$5.44\times {10}^{-16}$
$G_2$	1257	0.78	1.00	$3.85\times {10}^{-16}$	18610	32.58	0.94	$3.85\times {10}^{-16}$
$G_3$	674	0.42	1.00	$3.85\times {10}^{-16}$	731	1.61	1.00	$7.69\times {10}^{-16}$
$G_4$	8901	2.23	0.18	$1.48\times {10}^{-16}$	2317	4.38	0.22	$2.96\times {10}^{-16}$

Table 2. Table 2: H-eigenvalues of adjacency tensors of loose paths.

$m$	$r$	CSRH				CEST
		Iter.	Time(s)	Accu.	Err.	Iter.	Time(s)	Accu.	Err.
3	4	38123	9.14	1.00	$3.49\times {10}^{-16}$	42760	70.28	1.00	$3.49\times {10}^{-16}$
	6	62780	17.55	0.97	$5.67\times {10}^{-16}$	65706	105.53	0.99	$7.56\times {10}^{-16}$
	8	71311	23.38	0.66	$3.94\times {10}^{-16}$	76778	106.95	0.65	$7.88\times {10}^{-16}$
4	4	69517	16.92	1.00	$5.06\times {10}^{-16}$	49331	79.81	1.00	$5.06\times {10}^{-16}$
	6	86171	24.83	0.96	$5.55\times {10}^{-16}$	76105	113.11	0.98	$5.55\times {10}^{-16}$
	8	75907	24.71	0.33	$7.74\times {10}^{-16}$	91690	106.57	0.42	$9.68\times {10}^{-16}$

Table 3. Table 3: The p 𝑝 p -spectral radius of r 𝑟 r -uniform β 𝛽 \beta -stars.

n $p=3,r=3\left(p>r-1\right)$ Iter. Time(s) Accu. Err. 21 1835 0.34 1.00 $5.38\times {10}^{-16}$ 201 2609 0.60 1.00 $3.55\times {10}^{-15}$ 2,001 3539 1.87 1.00 $4.33\times {10}^{-14}$ 20,001 4475 12.93 1.00 $6.39\times {10}^{-14}$ 200,001 6038 263.39 0.98 $1.93\times {10}^{-11}$ 2,000,001 20018 15437.99 1.00 $1.22\times {10}^{-10}$ The $3$-spectral radius of $3$-uniform $\beta$-stars ( $p>r-1$ )

n Iter. Time(s) Accu. Err. 51 14747 4.79 0.99 $1.59\times {10}^{-11}$ 501 26019 14.52 0.98 $9.56\times {10}^{-12}$ 5,001 30108 57.82 0.99 $2.01\times {10}^{-11}$ 50,001 32387 426.60 0.95 $1.08\times {10}^{-11}$ 500,001 30070 6309.58 0.99 $4.49\times {10}^{-11}$ 5,000,001 51609 125869.02 0.97 $2.40\times {10}^{-10}$ The $4$-spectral radius of $6$-uniform $\beta$-stars ( )

Table 4. Table 4: p ϑ subscript 𝑝 italic-ϑ p_{\vartheta} -spectral radius of 3-uniform β 𝛽 \beta -star with 10 edges.

$p_n$	Iter.	Time(s)	Accu.	Err.
$p_{\vartheta }=\frac{12}{7}$	3037	0.99	1.00	$0.00$
$p_{\vartheta }=\frac{14}{9}$	13271	17.88	1.00	$3.08\times {10}^{-16}$
$p_{\vartheta }=\frac{10}{7}$	51018	110.53	1.00	$1.85\times {10}^{-16}$
$p_{\vartheta }=\frac{4}{3}$	84848	88.85	1.00	$3.07\times {10}^{-14}$

Table 5. Table 5: Top ten vertices in Figure 6 .

Ranking	$p=\frac{4}{3}$		$p=5$		$p=16$
	Num.	Val.	Num.	Val.	Num.	Val.
1	39	0.4082483175	41	0.4081204985	1	0.1709715830
2	38	0.4082482858	39	0.4081204985	31	0.1678396311
3	31	0.4082482855	31	0.4081204983	26	0.1618288319
4	41	0.4082482854	38	0.4081204982	39	0.1600192388
5	40	0.4082482849	40	0.4081204973	38	0.1600192387
6	37	0.4082482834	37	0.4081204958	41	0.1600192387
7	24	0.0000000000	28	0.0073198868	40	0.1600192386
8	34	0.0000000000	30	0.0073192175	37	0.1600192385
9	23	0.0000000000	26	0.0073061265	23	0.1550865094
10	3	0.0000000000	29	0.0071906282	22	0.1550865094

Table 6. Table 6: Top 10 authors.

Ranking	Author Name
	$p=2$	$p=12$	MultiRank
1	Zheng Chen	Wei-Ying Ma	C. Lee Giles
2	Wei-Ying Ma	Zheng Chen	Philip S. Yu
3	Qiang Yang	Jiawei Han	Wei-Ying Ma
4	Jun Yan	Philip S. Yu	Zheng Chen
5	Benyu Zhang	C. Lee Giles	Jiawei Han
6	Hua-Jun Zeng	Jian Pei	Christos Faloutsos
7	Weiguo Fan	Christos Faloutsos	Bing Liu
8	Wensi Xi	Yong Yu	Johannes Gehrke
9	Dou Shen	Qiang Yang	Gerhard Weikum
10	Shuicheng Yan	Ravi Kumar	Elke A. Rundensteiner

Equations

$${\mathbb{R}}^{[r,n]}\equiv {\mathbb{R}}^{\stackrel{r\text{-times}}{\overbrace{n\times n\times \cdots \times n}}}.$$

$$\mathcal{T}{\mathbf{x}}^r\equiv \sum _{i_1=1}^n\cdots \sum _{i_r=1}^n{t}_{i_1\cdots i_r}{\mathbf{x}}_{i_1}\cdots {\mathbf{x}}_{i_r}$$

$$(\mathcal{T}{\mathbf{x}}^{r-1}{)}_i\equiv \sum _{i_2=1}^n\cdots \sum _{i_r=1}^n{t}_{i{i}_2\cdots i_r}{\mathbf{x}}_{i_2}\cdots {\mathbf{x}}_{i_r},\text{for}i=1,\dots ,n.$$

$$\mathcal{T}{\mathbf{x}}^{m-1}=\lambda {\mathbf{x}}^{[m-1]},$$

$$\{\begin{array}{rlr}\mathcal{T}{\mathbf{x}}^{m-1}& =& \lambda \mathbf{x}\\ {\mathbf{x}}^{\top }\mathbf{x}& =& 1,\end{array}$$

$$w(G,\mathbf{x})=\sum _{e=\{i_1,\dots ,i_r\}\in E}s(e){\mathbf{x}}_{i_1}\cdots {\mathbf{x}}_{i_r},$$

$${\lambda }^{(p)}(G)=r!\underset{\Vert \mathbf{x}{\Vert }_p=1}{\max }w(G,\mathbf{x}),$$

$$\lambda_{L}(G)=\left\{\begin{array}[]{ll}\max&\hbox{$w(G,{\bf x})$}\\ \mathrm{s.t.}&\hbox{$\sum_{i=1}^{r}{\bf x}_{i}=1,$}\\ &\hbox{${\bf x}_{i}\geq 0,\quad\text{for}\,\quad i=1,\ldots,r.$}\end{array}\right.$$

$$a_{i_1\cdots i_r}=\{\begin{array}{rlrl}& \frac{s(e)}{(r-1)!}& & \text{if }\{i_1,\dots ,i_r\}\in E,\\ & 0& & \text{ otherwise. }\end{array}$$

$$\lambda_{L}(G)=\biggl{(}\begin{array}[]{c}n\\ r\end{array}\biggr{)}\frac{1}{n^{r}}.$$

$${\lambda }^{(p)}(G)=\underset{\Vert \mathbf{x}{\Vert }_p=1}{\max }(r-1)!\mathcal{A}{\mathbf{x}}^r$$

$${\lambda }^{(p)}(G)=\underset{\mathbf{x}\ne 0}{\max }(r-1)!\frac{\mathcal{A}{\mathbf{x}}^r}{\Vert \mathbf{x}{\Vert }_p^r}.$$

$$\{\begin{array}{rl}\max f(\mathbf{x})& =(r-1)!\frac{\mathcal{A}{\mathbf{x}}^r}{\Vert \mathbf{x}{\Vert }_p^r}\\ \text{s.t.}\Vert \mathbf{x}{\Vert }_2& =1.\end{array}$$

$$\nabla f(\mathbf{x})=\frac{r!}{\Vert \mathbf{x}{\Vert }_p^r}\left(\mathcal{A}{\mathbf{x}}^{r-1}-\mathcal{A}{\mathbf{x}}^r\Vert \mathbf{x}{\Vert }_p^{-p}{\mathbf{x}}^{⟨p-1⟩}\right),$$

$${\mathbf{x}}^{\top }\nabla f(\mathbf{x})=0$$

$$\underset{\vartheta \to \infty }{\lim }{p}_{\vartheta }=p_{\ast },$$

$$\underset{\vartheta \to \infty }{\lim }{\lambda }^{(p_{\vartheta })}(G)={\lambda }^{(p_{\ast })}(G).$$

$$\hat{f}(\mathbf{x},p)=(r-1)!\frac{\mathcal{A}{\mathbf{x}}^r}{\Vert \mathbf{x}{\Vert }_p^r}(\mathbf{x},p)\in {\mathbb{S}}^{n-1}\times (0,+\infty ),$$

$${\lambda }^{(p)}(G)=\underset{\mathbf{x}\in {\mathbb{S}}^{n-1}}{\max }\hat{f}(\mathbf{x},p).$$

$$\hat{f}({\mathbf{x}}_{\vartheta }^{\ast },p_{\vartheta })={\lambda }^{(p_{\vartheta })}(G).$$

$$\underset{\vartheta \to \infty }{\lim }{\mathbf{x}}_{\vartheta }^{\ast }={\mathbf{x}}_0^{\ast }.$$

$$\hat{f}(\tilde{\mathbf{x}},p_{\vartheta })\le \hat{f}({\mathbf{x}}_{\vartheta }^{\ast },p_{\vartheta })$$

$$\underset{\vartheta \to \infty }{\lim }\hat{f}(\tilde{\mathbf{x}},p_{\vartheta })\le \underset{\vartheta \to \infty }{\lim }\hat{f}({\mathbf{x}}_{\vartheta }^{\ast },p_{\vartheta }).$$

$$\hat{f}(\tilde{\mathbf{x}},p_{\ast })\le \hat{f}({\mathbf{x}}_0^{\ast },p_{\ast })$$

$$\hat{f}({\mathbf{x}}_0^{\ast },p_{\ast })=\underset{\vartheta \to \infty }{\lim }\hat{f}({\mathbf{x}}_{\vartheta }^{\ast },p_{\vartheta })=\underset{\vartheta \to \infty }{\lim }{\lambda }^{(p_{\vartheta })}(G),$$

$${\mathbf{d}}_k^{\top }\nabla f({\mathbf{x}}_k)>0.$$

$$({\mathbf{x}}_{k+1}+{\mathbf{x}}_k{)}^{\top }{\mathbf{d}}_k=0.$$

$$({\mathbf{x}}_k+{\mathbf{x}}_{k+1}{)}^{\top }{W}_k({\mathbf{x}}_{k+1}+{\mathbf{x}}_k)=-({\mathbf{x}}_{k+1}+{\mathbf{x}}_k{)}^{\top }{W}_k({\mathbf{x}}_{k+1}+{\mathbf{x}}_k)=0.$$

$${\mathbf{d}}_k=W_k({\mathbf{x}}_k+{\mathbf{x}}_{k+1}).$$

$${\mathbf{p}}_k^{\top }\nabla f({\mathbf{x}}_k)>0.$$

Full text

Key words.   Eigenvalue, hypergraph, large scale tensor, network analysis, pagerank, $p$-spectral radius.

AMS subject classifications.   05C65, 15A18, 15A69, 65F15, 65K05, 90C35, 90C53

1 Introduction

With the emergence of big data in various field of our social life, it becomes significant and challenging to analyze the massive data and extract valuable information from them. Hypergraph, as an extension of graph, provides an efficient way to represent complex relationships among objects in applied science, such as chemistry [37, 33], computer science [24, 56, 30], and image processing [5, 19, 11]. The spectral hypergraph theory has been widely studied in [14, 29, 32, 42, 53, 65, 67], which reveal combinatorial and geometric structures of hypergraphs. Moreover, spectral hypergraph approaches are useful tools to address issues in real world. Spectral hypergraph partitioning and spectral hypergraph clustering have broad applications in network analysis [43, 59], image segmentation [18], multi-label classification [61], machine learning [68], and data analysis [2, 40]. Hypergraph spectral hashing techniques highly contribute to problems of similarity search and retrieval of social image [69, 41].

In this paper, we focus on the computation of $p$-spectral radii of uniform hypergraphs. The $p$-spectral radius of a hypergraph was introduced in [32] and linked with extremal hypergraph problems. Extremal graph theory, as a branch of graph theory, is one of the most attractive and best studied area in combinatorics. Turán [63] introduced the famous Turán graph and Turán theorem in 1941, when is regarded as the start of the extremal graph theory. Naturally, the question was extended from graph to hypergraph in [64] that is to find the largest number of edges in a hypergraph which is $F$-free111 A uniform hypergraph that does not have a subgraph isomorphic to the uniform hypergraph $F$ is said to be $F$-free.. Although the Turán-type problem is adequately complete for ordinary graphs, cases are much more challenging when it comes to hypergraph. In [48], Nikiforov proved the spectral Turán-type inequality which generalized the Turán theorem. In [32], the $p$-spectral version of Nikiforov’s inequality and the $p$-spectral version of a hypergraph Turán result were given, and it was showed that this result can be employed in solving ‘degenerate’ Turán-type problems. Furthermore, it was proved that the edge extremal problems are asymptotically equivalent to the extremal $p$-spectral radius problems in [49].

The $p$-spectral radius of a hypergraph covers not only the number of edges in extremal problems, but also the notions, such as Lagrangian, and the spectral radius of a hypergraph [42]. When $p=1$, the $p$-spectral radius of a hypergraph turns out to be its Lagrangian. The Lagrangians of graph and hypergraph were proposed in [44] to prove the Turán’s theorem for graphs. The Largrangians of hypergraphs were used to disprove the conjecture of Erdös [20, 22] and to find non-jumping numbers for hypergraphs [23, 54, 55]. Also, the Lagrangian of a hypergraph is associated with problems of determining Turán densities of hypergraphs [6, 31, 45, 60], which is an asymptotic solution to a (non-degenerate) Turán problem. When $p=2$, the $p$-spectral radius of a uniform hypergraph is the largest Z-eigenvalue [57] of its adjacency tensor. When $p$ is even and equals the order of this hypergraph, the $p$-spectral radius becomes the largest H-eigenvalue of the adjacency tensor of $G$. Therefore, the $p$-spectral radius is connected with the (adjacency) spectral radius of a hypergraph [27, 39, 42]. Additionally, Kang et al. provided solutions to several $p$-spectral radius related extremal problems in [29]. Nikiforov in [50] did a comprehensive study and obtained many theoretical conclusions about $p$-spectral radius.

Apart from the application in extremal hypergraph theory, the $p$-spectral radius model constructs a framework to quantify the importance of objects or centrality in networks. Evaluating the significance or popularity of objects is a significant problem in data mining. It can be used to determine the importance of web pages [52, 36, 16], forecast customer behaviour [38], retrieve images [28] and so on. In the $p$-spectral radius model, entries of the vector associated with the $p$-spectral radius of a hypergraph are called $p$-optimal weighting and represent the significance of its corresponding vertices. The ranking result varies when $p$ changes. We will explain the meaning of different ranking results and show the numerical performance of our algorithm in sorting real-life data in Section 6.

Calculation of $p$-spectral radii of hypergraphs is related to several methods for evaluating tensor eigenvalues. Algorithms for tensor eigenvalues, such as the shifted symmetric higher-order power method (SS-HOPM) in [34], the generalized eigenproblem adaptive power (GEAP) method in [35], an extension of Collatz’s method (NQZ) in [46], and the CEST method, can be employed when $p$ equals 2 or when $p$ equals the order of an even-uniform hypergraph. When $p$ is even, the $p$-spectral radius problem is equivalent to the generalized tensor eigenvalue problem [9, 17]. Therefore, methods for this generalized tensor eigenvalue problem, such as the polynomial optimization related algorithm for finding all real eigenvalues of a symmetric tensor given by Cui et al. in [15], and the homotopy approach for all eigenpairs of general real or complex tensors proposed by Chen et al. in [10] can be employed to compute even $p$-spectral radius of small scale hypergraphs. However, the problem of computing $p$-spectral radii of arbitrary hypergraph is still open. This is the main motivation of our paper.

To solve the $p$-spectral radius problem, we introduce a spherically constrained maximization model, which is equivalent to the original problem. Then we use an effective conjugate gradient method to acquire an ascent direction for the constrained optimization model. Next, we employ the Cayley transform to project the ascent direction on the unit sphere. It is proved that there exists a positive parameter in the curvilinear line search such that the Wolfe conditions hold. Based on the above foundation, we propose a numerical method for computing $p$-spectral radii of hypergraphs (CSRH) with $p>1.$ When $p=1,$ the CSRH method is able to approximate the $1$-spectral radii (Largrangians) of hypergraphs. In the convergence analysis, we prove that the CSRH algorithm is convergent and it converges to the global optimization point with high probability. Numerical experiments show that CSRH is preponderant when compared to existing methods for computing Z-eigenvalues and H-eigenvalues of adjacency tensors. Moreover, CSRH is capable of calculating $p$-spectral radii of hypergraphs with millions of vertices effectively. In addition, we find that the significance of vertices of hypergraphs is related to the order of elements of the $p$-optimal weighting. Therefore, we apply the CSRH method to rank the vertices of the corresponding hypergraph from different viewpoints when $p$ is different, which is useful in network analysis. As an example, we show that our numerical results agree with the observed data of a small weighted hypergraph. Furthermore, we successfully rank 10305 authors based on their publication information by establishing a hypergraph model and using CSRH to solve the corresponding $p$-spectral radius problem. We sort the authors from the view of individual and group respectively. The result of our ranking can be reasonably explained and are in line with the existing consequences in [47].

The paper is organized as follows. In Section $2$, we introduce mathematical notions. The computational issues about $p$-spectral radius are addressed in Section $3$, where our new method CSRH for computing $p$-spectral radii of hypergraphs is given. In Section $4$, we analyze the convergent property of the CSRH method. The numerical experiments are represented in Section $5.$ In Section $6,$ we show the application of CSRH method in network analysis. The ranking results of a toy example and a large scale real-world problem are presented. Finally, we draw conclusions in Section $7$.

2 Preliminary

In this section we introduce useful notions and important results on hypergraphs and tensors. Let $\mathbb{R}^{[r,n]}$ be the $r$th order $n$-dimensional real-valued tensor space, i.e.,

[TABLE]

A tensor $\mathcal{T}=(t_{i_{1}\cdots i_{r}})\in\mathbb{R}^{[r,n]}$ with $i_{j}=1,\ldots,n$ for $j=1,\ldots,r,$ is said to be symmetric, if $t_{i_{1}\cdots i_{r}}$ is unchanged under any permutation of indices [13]. Two operations between $\mathcal{T}$ and any vector ${\bf x}\in\mathbb{R}^{n}$ are defined as

[TABLE]

and

[TABLE]

Note that, $\mathcal{T}{\bf x}^{r}\in\mathbb{R}$ and $\mathcal{T}{\bf x}^{r-1}\in\mathbb{R}^{n}$ are a scalar and a vector respectively, and $\mathcal{T}{\bf x}^{r}={\bf x}^{\top}(\mathcal{T}{\bf x}^{r-1}).$

If there exists a real number $\lambda$ and a nonzero real vector ${\bf x}$ satisfying

[TABLE]

then $\lambda$ is called an H-eigenvalue of $\mathcal{T}$ with ${\bf x}$ being the associated H-eigenvector [57, 58]. Additionally, ${\bf x}^{[m-1]}\in\mathbb{R}^{n}$ is a vector, of which the $i$th element is ${\bf x}_{i}^{m-1}.$ When a real vector ${\bf x}$ and a real number $\lambda$ satisfy the following system

[TABLE]

$\lambda$ is called a Z-eigenvalue of $\mathcal{T}$ and ${\bf x}$ is the corresponding Z-eigenvector [57].

Definition 2.1** (Hypergraph).**

A hypergraph is defined as $G=(V,E)$, where $V=\{1,2,\ldots,n\}$ is the vertex set and $E=\{e_{1},e_{2},\ldots,e_{m}\}\subseteq 2^{V}$ (the powerset of $V$) is the edge set. We call $G$ an $r$-uniform hypergraph when $|e_{p}|=r\geq 2$ for $p=1,\ldots,m$ and $e_{i}\neq e_{j}$ in case of $i\neq j.$

If each edge of a hypergraph is linked with a positive number $s(e),$ then this hyperpragh is called a weighted hypergraph and $s(e)$ is the weight associated with the edge $e.$ An ordinary hypergraph can be regarded as a weighted hypergraph with the weight of each edge being $1.$

In the rest of this paper, an $r$-uniform hypergraph is abbreviated to an $r$-graph for convenience and hence the hypergraph $G$ refers to an $r$-graph. The degree of a vertex $i\in V$ is given by $d(i)=\mathrm{sum}\{s(e):i\in e,e\in E\}.$ The weight polynomial of $G$ [62] is defined as

[TABLE]

in which ${\bf x}$ is a vector in $\mathbb{R}^{n}$, $e=\{i_{1},\ldots,i_{r}\}$ is an edge of $G$ and $s(e)$ is the weight of $e.$

Definition 2.2** ($p$-spectral radius [32, 29]).**

When $p\geq 1$, the $p$-spectral radius of $G$, denoted by $\lambda^{(p)}(G)$, is defined as

[TABLE]

and we call any vector ${\bf x}$ solving (2.3) a $p$-optimal weighting of $G$ [7].

When $p=1,$ the $p$-spectral radius of $G$ coincides with its Lagrangian $\lambda_{L}(G)$ [21, 62], which is defined as

[TABLE]

The vector ${\bf x}$ related to the Lagrangian of $G$ is named the optimal legal weighting [7, 62].

Definition 2.3** (Adjacency tensor ).**

The adjacency tensor $\mathcal{A}$ of a weighted $r$-graph $G$ is defined as an $r$th order $n$-dimensional symmetric tensor with its elements being

[TABLE]

It is obvious from (2.3) that the $2$-spectral radius is exactly the product of $(r-1)!$ times the largest Z-eigenvalue of the adjacency tensor $\mathcal{A}$, and when $r$ is even the $r$-spectral radius is $(r-1)!$ times the largest H-eigenvalue of $\mathcal{A}$ [57].

Although there is no general formula or algorithm for us to compute the $p$-spectral radius of a hypergraph directly, research on $p$-spectral radius of hypergraphs with certain structures has made some progress.

Theorem 2.1** ([50]).**

*Let $r$-graph $G$ be a $\beta$-star with $m$ edges .

a. If $p>r-1,$ then $\lambda^{(p)}(G)=r!r^{-\frac{r}{p}}m^{(1-\frac{r-1}{p})}.$

b. If $p<r-1,$ then $\lambda^{(p)}(G)=r!r^{-\frac{r}{p}}.$

c. If $p=r-1,$ then $\lambda^{(p)}(G)=(r-1)!r^{-\frac{1}{r-1}}.$*

Proposition 2.1** ([7]).**

If $G$ is a complete $r$-graph with $n$ vertices, then the Lagrangian of $G$ is

[TABLE]

A multiset is an extension of the ordinary set, such that the objects or elements in the multiset are repeatable. If the edge set $E$ of a hypergraph $G$ is a set of multisets, then $G$ is called a multi-hypergraph [53]. Naturally, the $p$-spectral radius problem can be extended from hypergraph to muli-hypergraph. The algorithm and theoretical analysis in the following part of this paper are also applicable to $p$-spectral radius problems of multi-hypergraphs. In the rest of this paper, the symbol $\|\cdot\|$ refers to $\ell_{2}$ norm and the parameter $p$ is a positive integer unless stated otherwise.

3 Computation of the $p$-spectral radius of a hypergraph

We transform the $p$-spectral radius in (2.3) into a spherically constraint optimization problem and propose an iterative algorithm to solve it.

3.1 Spherically constraint form for $\lambda^{(p)}(G)$

The $p$-spectral radius of $G$ in (2.3) can be reformulated as

[TABLE]

where $\mathcal{A}$ is the adjacency tensor of $G.$ The maximization problem (3.1) is equivalent to an unconstrained format, that is

[TABLE]

In order to restrict the search region and keep the vector ${\bf x}$ away from zero, we add a spherically constraint on $\lambda^{(p)}(G)$ in (3.2). Due to the zero-order homogeneous property of $\mathcal{A}{\bf x}^{r}/\|{\bf x}\|_{p}^{r},$ we can obtain $\lambda^{(p)}(G)$ by solving the following problem

[TABLE]

When $p>1$, the objective function $f({\bf x})$ is differentiable for any nonzero ${\bf x}$ and the gradient of $f({\bf x})$ is

[TABLE]

where ${\bf x}^{\langle p-1\rangle}$ represents a vector whose $i$th element is $({\bf x}^{\langle p-1\rangle})_{i}=|x_{i}|^{p-1}\text{sgn}(x_{i}).$ Since $f({\bf x})$ is zero-order homogeneous, we have

[TABLE]

for any $0\neq{\bf x}\in\mathbb{R}^{n}.$

Based on the spherically constrained form in (3.3), we have the following proposition, which provides a way to approximate the $p$-spectral radius of a hypergraph when it cannot be computed directly.

Proposition 3.1**.**

Let $p_{\vartheta}$ be a sequence such that

[TABLE]

where each $p_{\vartheta}>0.$ Then

[TABLE]

Proof.

We restrict the domain of ${\bf x}$ on a unit sphere, which is denoted as $\mathbb{S}^{n-1}\equiv\{{\bf x}\in\mathbb{R}^{n}:{\bf x}^{\top}{\bf x}=1\}.$ Rename the function in (3.3) as

[TABLE]

and we have

[TABLE]

Here $\hat{f}({\bf x},p)$ is continuous. Let $\{{\bf x}_{\vartheta}^{*}\}$ be an infinite sequence on the compact space $\mathbb{S}^{n-1},$ such that

[TABLE]

If there are more than one point satisfying the equation (3.8), we randomly choose one of them to be ${\bf x}_{\vartheta}^{*}.$ Suppose $\{{\bf x}_{\vartheta}^{*}\}$ is a convergent sequence without loss of generality. Since the sequence is bounded, there exists a point ${\bf x}_{0}^{*}\in\mathbb{S}^{n-1}$ satisfying

[TABLE]

For any $\tilde{{\bf x}}\in\mathbb{S}^{n-1},$ we have

[TABLE]

from (3.8), which indicates that

[TABLE]

Then we obtain

[TABLE]

based on (3.6) and (3.9). Therefore we have $\hat{f}({\bf x}_{0}^{*},p_{*})=\max_{{\bf x}\in\mathbb{S}^{n-1}}\hat{f}({\bf x},p_{*})=\lambda^{(p_{*})}(G).$ Since

[TABLE]

conclusion (3.7) is then obtained.

∎

3.2 The CSRH algorithm

We employ an iterative algorithm to solve (3.3).

Suppose that the current iterate is a unit vector ${\bf x}_{k}$. Our task is to find a new iterate ${\bf x}_{k+1},$ which satisfies the following two conditions.

${\bf x}_{k+1}$ is on the unit sphere; 2. 2.

${\bf d}_{k}={\bf x}_{k+1}-{\bf x}_{k}$ is an ascent direction, i.e.,

[TABLE]

In Figure 1, the current iterate ${\bf x}_{k}$ is on the unit sphere and we can see that ${\bf x}_{k+1}$ is a unit vector if and only if the vector ${\bf x}_{k+1}+{\bf x}_{k}$ and the vector ${\bf d}_{k}={\bf x}_{k+1}-{\bf x}_{k}$ are perpendicular to each other, i.e.

[TABLE]

Let $W_{k}$ be a skew-symmetric matrix, i.e., $W_{k}=-W_{k}^{\top}.$ Then we have

[TABLE]

Therefore, the equation (3.13) is feasible and the first condition of ${\bf x}_{k+1}$ holds when

[TABLE]

Furthermore, based on the optimization techniques it is available to find an ascent direction ${\bf p}_{k}$ such that

[TABLE]

Then the existing information in Figure 1 for us to obtain ${\bf d}_{k}$ is ${\bf p}_{k}$ and ${\bf x}_{k},$ both of which have relation with $\nabla f({\bf x}_{k})$ in (3.15) and (3.5) respectively. Hence, in order to satisfy (3.12) we construct ${\bf d}_{k}$ as a combination of ${\bf x}_{k}$ and ${\bf p}_{k},$ i.e.,

[TABLE]

and obtain

[TABLE]

Therefore, if $b>0$ in (3.16), ${\bf d}_{k}$ is an ascent direction with ${\bf d}_{k}^{\top}\nabla f({\bf x}_{k})>0$.

The previous analysis shows that the two conditions of ${\bf x}_{k+1}$ are valid when ${\bf d}_{k}$ satisfies (3.14) and (3.16) for $b>0$. This motivates us to construct the skew-symmetric matrix $W_{k}$ by ${\bf x}_{k}$ and ${\bf p}_{k}.$ Let

[TABLE]

with $\alpha$ being a positive parameter. The constant $b=\frac{1}{2}\alpha{\bf x}_{k}^{\top}({\bf x}_{k}+{\bf x}_{k+1})$ in (3.16). Since the angle between vectors ${\bf x}_{k}$ and ${\bf x}_{k}+{\bf x}_{k+1}$ is less than or equal to $\frac{\pi}{2}$ in Figure 1, then we have $b\geq 0.$ However if $b=0,$ i.e., ${\bf x}_{k+1}=-{\bf x}_{k},$ there is a contradiction when we substitute ${\bf x}_{k+1}$ by $-{\bf x}_{k}$ in (3.14). Hence, we have $b>0$ and equations (3.14) and (3.16) hold, which means the two conditions of ${\bf x}_{k+1}$ are satisfied when $W_{k}$ is the matrix in (3.18) with ${\bf p}_{k}$ being an ascent direction.

Lemma 3.1**.**

The new iterate ${\bf x}_{k+1}$ can be expressed as

[TABLE]

from (3.14) and (3.18). Further we have

[TABLE]

Proof.

From (3.14), we obtain ${\bf x}_{k+1}=Q{\bf x}_{k},$ where

[TABLE]

That is to say the orthogonal transform is in fact the Cayley transform. The proof is then similar to Lemma $3.2$ in [8, 12]. ∎

For the new point ${\bf x}_{k+1}$ in (3.19), a crucial step is to find an ascent direction ${\bf p}_{k}$ to guarantee the ascent property in (3.15). Since problems related with hypergraphs and tensors are often large and time-consuming for computation, we employ the nonlinear conjugate gradient method, which is proposed for large-scale nonlinear optimization problems, to acquire a suitable ${\bf p}_{k}$. The nonlinear conjugate gradient method does not need the Hessian matrices of the objective function and is usually faster than the steepest descent method. In [25, 26], a nonlinear conjugate gradient method called CG$\_$DESCENT was given and it was proved that the CG$\_$DESCENT possesses a good descent property. Attracted by this merit, we adopt the construction of parameter $\beta_{k}$ in CG$\_$DESCENT and obtain the ascent direction ${\bf p}_{k}$ by

[TABLE]

The scalar $\beta_{k-1}$ above is defined as $\beta_{k-1}=\max(0,\tilde{\beta}_{k-1})$, where

[TABLE]

$\textbf{y}_{k-1}=\nabla f({\bf x}_{k})-\nabla f({\bf x}_{k-1}),$ parameters $\frac{1}{4}<\tau<1$ and $\epsilon>0.$ The initial direction is chosen as ${\bf p}_{0}=\nabla f({\bf x}_{0}).$ The direction ${\bf p}_{k}$ in (3.21) is proved to satisfy the ascent property in the following Lemma.

Lemma 3.2**.**

The search direction ${\bf p}_{k}$ generated by (3.21) satisfies the sufficient ascent condition, i.e.

[TABLE]

and there exists a constant $M_{0}>1$ such that

[TABLE]

Proof.

When $\beta_{k}=0,$ it is easy to show that the two inequalities hold. For $\beta_{k}\neq 0,$ we have

[TABLE]

Since

[TABLE]

we obtain

[TABLE]

Then we deduce that

[TABLE]

Inequality (3.24) is valid when $M_{0}=1+\frac{1}{\epsilon}+\frac{\tau}{\epsilon^{2}}.$ ∎

In the curvilinear line search, the parameter $\alpha$ in (3.19) is determined to ensure that the Wolfe conditions hold. We provide the details in the next subsection.

3.3 Feasibility of Wolfe conditions

In this section we prove that there exists a step length $\alpha_{k}$ satisfying the Wolfe conditions for the curvilinear search in (3.19) in each iteration. First, we compute the derivative of $\alpha$ which plays an important role in line search.

Lemma 3.3**.**

Let $f^{\prime}(\alpha)$ be the derivative of $f({\bf x}_{k+1}(\alpha))$ at point $\alpha.$ Then we have

[TABLE]

Proof.

Equation (3.19) means that

[TABLE]

Then we take derivative with respect to $\alpha$ as follows

[TABLE]

By multiplying both sides of (3.26) by $\alpha$ we get

[TABLE]

from (3.19). Since $\nabla f({\bf x}_{k+1}(\alpha))^{\top}{\bf x}_{k+1}(\alpha)=0,$ from (3.27) we obtain

[TABLE]

∎

Since $f({\bf x})$ is twice continuously differentiable in the compact set $\mathbb{S}^{n-1}$ , we can find a constant $M$ such that

[TABLE]

For a given optimization algorithm which enjoys a good ascent or descent property, it is proved that step lengths that satisfy the Wolfe conditions exist for a monotonous line search in [51, Lemma 3.1]. In the following theorem we prove that Wolfe conditions are practicable for the curvilinear line search in our algorithm.

Theorem 3.1**.**

If  $0<c_{1}<c_{2}<1$, there exists $\alpha_{k}>0$ satisfying

[TABLE]

Proof.

Let ${\bf x}(\alpha)={\bf x}_{k+1}(\alpha)$ and $f(\alpha)=f({\bf x}_{k+1}(\alpha)).$ From (3.19), we have ${\bf x}_{k+1}^{\prime}(0)=-{\bf x}_{k}^{\top}{\bf p}_{k}{\bf x}_{k}+{\bf p}_{k},$ and

[TABLE]

Denote a linear function $l(\alpha)=f({\bf x}_{k})+c_{1}\alpha\nabla f({\bf x}_{k})^{\top}{\bf p}_{k}.$ Then $f(0)=l(0)=f({\bf x}_{k})$ and $f^{\prime}(0)>l^{\prime}(0)>0$ due to $0<c_{1}<1$ and $\nabla f({\bf x}_{k})^{\top}{\bf p}_{k}>0$ in (3.23). Since $f(\alpha)$ is bounded above, the graph of $f(\alpha)$ must intersect with the line $l(\alpha)$ at least once when $\alpha>0$. Suppose $\bar{\alpha}$ is the smallest intersection point, we obtain

[TABLE]

By the mean value theorem, we can find $\rho\in(0,\bar{\alpha})$ satisfying

[TABLE]

On the other hand, from (3.5) and (3.19) we have

[TABLE]

Then we have

[TABLE]

Combining (3.32) and (3.33), we have

[TABLE]

Further, from (3.31) we obtain

[TABLE]

Combing (3.34) and (3.35) we have

[TABLE]

Since

[TABLE]

and $|{\bf x}_{k}^{\top}{\bf p}_{k}|\leq\|{\bf p}_{k}\|$, we have

[TABLE]

Since $\nabla f({\bf x}_{k})^{\top}{\bf p}_{k}\geq 0$,

[TABLE]

Since $c_{2}>c_{1}$, inequality (3.30) holds when $\alpha_{k}=\rho.$ Also from the condition $\rho\in(0,\bar{\alpha})$, we have $f(\alpha_{k})>l(\alpha_{k})$ and (3.29) is obtained. ∎

Up to now, the algorithm CSRH for computing the $p$-spectral radius of a hypergraph is available. First we transform the original model of $\lambda^{(p)}(G)$ into an equivalent constrained optimization problem on the unit sphere (3.3). To solve the constrained model, we compute the ascent direction ${\bf p}_{k}$ from (3.4), (3.22) and (3.21), and choose a proper $\alpha_{k}$ so that the next iterate gained via (3.19) satisfies the Wolfe conditions (3.29) and (3.30). A fast computation method for calculating $\mathcal{A}{\bf x}^{r}$ and $\mathcal{A}{\bf x}^{r-1}$ was proposed in [8], which improves the efficiency of products of adjacency tensor and vector. We also adopt this technique in our algorithm.

4 Convergence analysis

In this section we prove that the CSRH algorithm converges to a stationary point of $f({\bf x})$ and touches the exact $p$-spectral radius with a high probability. Our CSRH algorithm terminates finitely when there exits a constant $c$ such that $\nabla f({\bf x}_{c})=0.$ The following convergence analysis is for the case that the sequence $\{{\bf x}_{k}\}$ is infinite and $\nabla f({\bf x}_{k})$ is always a nonzero vector.

4.1 Convergence results

Next theorem shows that CSRH algorithm is convergent.

Theorem 4.1**.**

Suppose the sequence $\{{\bf x}_{k}\}$ is generated by the algorithm CSRH from any ${\bf x}_{0}\in\mathbb{S}^{n}$. Then we have

[TABLE]

Proof.

The demonstration is divided into two steps. First, we show that the Zoutendijk condition holds, i.e.,

[TABLE]

Here $\varphi_{k}$ is the angle between $\nabla f({\bf x}_{k})$ and ${\bf p}_{k}$, which is denoted as

[TABLE]

Since $\nabla^{2}f({\bf x})$ is bounded, we have $\nabla f({\bf x})$ is Lipschitz continuous on $\mathbb{S}^{n-1}$, i.e.,

[TABLE]

for a constant $L>0$. From (3.18), we have

[TABLE]

Hence from (3.14)

[TABLE]

From (4.2) and (4.3), we have

[TABLE]

From (3.30), we obtain

[TABLE]

By using the above two relations, we can derive the inequality

[TABLE]

which implies

[TABLE]

Then from (3.29), we obtain

[TABLE]

which derives the following inequality

[TABLE]

Since $f({\bf x})$ is bounded in (3.28), the inequality (4.1) is then deduced.

Next, we show that the angle $\varphi_{k}$ is bounded away from $\frac{\pi}{2}$. By combining (3.23) and (3.24), we obtain

[TABLE]

The above inequalities indicate that

[TABLE]

Therefore, from (4.1) we have

[TABLE]

∎

Recall that the graph of a function $h({\bf x})$ is defined as

[TABLE]

For the function $f({\bf x})$ involved in our problem (3.3), we have

[TABLE]

where $p$ and $r$ are positive integers. Since $\text{Gr}\,f$ is a semialgebraic set, $f({\bf x})$ is a semialgebraic function and satisfies the Łojasiewicz inequality [1, 4, 66], which means that for a critical point ${\bf x}_{*}$ of $f({\bf x})$, there exist constants $\theta\in[0,1)$ and $C_{1}>0$, as well as $\mathscr{U}$ being a neighbourhood of ${\bf x}_{*}$ such that

[TABLE]

for ${\bf x}\in\mathscr{U}.$ The next theorem shows that if the sequence $\{{\bf x}_{k}\}$ is infinite, it has a unique accumulation point.

Theorem 4.2**.**

Assume the infinite sequence $\{{\bf x}_{k}\}$ is generated by the CSRH algorithm. Then it converges to a unique point ${\bf x}_{*},$ that is,

[TABLE]

and ${\bf x}_{*}$ is a first-order stationary point.

Proof.

From (4.5), (3.23) and (3.24) we have

[TABLE]

Moreover, from (3.29) and (3.23) we obtain

[TABLE]

We take no account of condition $\|\nabla f({\bf x}_{k})\|=0$ under which the algorithm terminates finitely. The above inequality indicates that

[TABLE]

Based on (3.24), (3.29) and (4.3), we have

[TABLE]

From (4.9) and (4.10), as well as the Łojasiewicz inequality (4.7), we have the conclusions hold based on [1, Theorem 3.2]. ∎

4.2 Probability of obtaining the exact $p$-spectral radius

Due to the feasibility of Łojasiewicz inequality in (4.7), we get the probability of the CSRH method touching the true $p$-spectral radius.

Proposition 4.1** ( ).**

Suppose CSRH algorithm is implemented from $N$ uniformly distributed initial points on $\mathbb{S}^{n-1}$ for $N$ times. We take the largest one among the results of these trails as the $p$-spectral radius of the relevant problem. The probability of getting the exact $p$-spectral radius is

[TABLE]

in which $\zeta$ is a constant satisfying $\zeta\in(0,1].$ If $N$ is large enough, the probability is high.

Proof.

This Proposition can be proved in the way similar to [8, Theorem 4.9]. We omit the details. ∎

5 Numerical experiments

In this section, we show the performance of CSRH for computing $p$-spectral radii of both small and large scale hypergraphs. We compare our method with several existing methods for computing eigenvalues of adjacency tensors. Examples of approximating the Lagrangian of a hypergraph are given in Subsection 2. All experiments are carried out by using MATLAB version R2015b and Tensor Toolbox version 2.6 [3]. The experiments in Subsections 5.1 and 5.2 are terminated when

[TABLE]

where $\lambda^{(p)}$ is our computed $p$-spectral radius and $\lambda^{(p)}_{*}(G)$ is the exact result obtained from theorems or conclusions in existing literature. The maximum iteration of CSRH is taken as 1000 for all algorithms except those performed by the MATLAB function in Tensor Toolbox. For each experiment in this section, we compute 100 times to obtain 100 estimated values $\lambda^{(p)}_{1},\ldots,\lambda^{(p)}_{100}$ and choose the largest one as our computational result of the $p$-spectral radius related with $G$. When $\lambda^{(p)}_{*}(G)$ is attainable, the accuracy rate of the CSRH algorithm is defined as

[TABLE]

Each number of iterations (Iter.) and computational time (Time) we reported in this section is the sum of corresponding quantities for all 100 executions of the experiment. The relative errors (Err.) between the numerical results and the exact solutions are provided.

5.1 Computation of $p$-spectral radii of hypergraphs

We compare the following three algorithms for computing eigenvalues of adjacency tensors associated with different hypergraphs:

• •

An adaptive shifted power method [34] SS-HOPM. This method can be invoked by eig_sshopm in Tensor Toolbox 2.6 for Z-eigenvalues of symmetric tensors.

• •

A first-order optimization algorithm CEST [8] which is proposed for eigenvalues of large scale sparse tensors involving even order hypergraphs.

• •

CSRH: the method proposed in Section 3.

Example 1 ($\mathbf{p=2}$). First, we compute the largest Z-eigenvalues of adjacency tensors of the following hypergraphs:

[TABLE]

The first hypergraph $G_{1}$ is given in [65] as Example 1, while the last three hypergraphs are Example $4,$ $7$ and $9$ in [53]. The hypergraph $G_{4}$ is actually a tetrahedron.

In Table 1, we demonstrate results of CSRH and SS-HOPM for computing the largest Z-eigenvalues of adjacency tensors of some small hypergraphs. Since all the four hypergraphs given above are of odd orders, the comparison does not include CEST method, which is designed for even order hypergraphs. The Err. column shows the relative error between the computational result and the exact largest Z-eigenvalue provided in the corresponding references. Under the condition that the relative error reaches $10^{-16}$, our CSRH method is much more stable and efficient than the SS-HOPM method.

In the next experiment, we study the probability of CSRH method getting the true largest Z-eigenvalue of $G_{4}$ and show that the probability increases along with the trail times. We employ the CSRH method to compute the largest Z-eigenvalue of the adjacency tensor of $G_{4}$ from uniformly distributed and randomly chosen initial points. Once the relative error between the computational largest Z-eigenvalue and its exact value $3/2$ reaches $10^{-8},$ the experiment is terminated and we record the number of trails. This experiment is repeated for one thousand times. Let $\sigma(i)$ be the total occurrence of experiments whose trail time is the integer $i.$ The frequency of touching the exact Z-eigenvalue when running $i$ times is

[TABLE]

In Figure 2, we display the relation between trail times and success probability. It illustrates that the probability tends to one along with the increase of trail times $i,$ which coincides with the conclusion in Theorem 4.1.

Example 2 ($\mathbf{p=r}$). Next, we compare CEST and CSRH methods for computing the largest H-eigenvalues of adjacency tensors of loose paths. An $r$-graph with $m$ edges is called a loose path if its vertex set is

[TABLE]

and its edge set is

[TABLE]

An $r$-uniform loose path with $m$ edges has $m(r-1)+1$ vertices. For example, the $6$-unform loose path with $4$ edges in Figure 3 has $21$ vertices.

The following theorem proved in [67] offers a convenient way to acquire the largest H-eigenvalues of adjacency tensors of loose paths with $m=3$ or $m=4$.

Theorem 5.1** ([67]).**

Let $G$ be an $r$-uniform loose path with $m$ edges and $\lambda_{H}(G)$ be the largest H-eigenvalue of its adjacency tensor $\mathcal{A}$. Then we have

\lambda_{H}(G)=\big{(}\frac{1+\sqrt{5}}{2}\big{)}^{\frac{2}{r}}* for $m=3,$* 2. 2.

$\lambda_{H}(G)=3^{\frac{1}{r}}$* for $m=4.$*

In Table 2, we compare CSRH and CEST for computing the largest H-eigenvalues of adjacency tensors of different loose paths. The column Err. presents the relative error between our computed result and the exact one given by Theorem 5.1. When relative error achieves precision of $10^{-16},$ the CSRH method saves at least $75\%$ of the time CEST takes in every problem. The comparison between CEST and CSRH verifies that the high efficiency of CSRH method does not only relies on the fast computation technique in [8], because CEST method use this technique as well.

Example 3. If all edges of a hypergraph share a same vertex, then it is called a $\beta$-star. An $r$-uniform $\beta$-star with $m$ edges have $m(r-1)+1$ vertices.

We present a class of $6$-uniform $\beta$-star in Figure 4 as an example.

We calculate $p$-spectral radii of $\beta$-stars with various orders and edges and display the results in Table 3. The Err. column presents the relative error between our computational result and the corresponding exact result generated from Theorem 2.1. It can be seen that all tests succeed with high accuracy rates. Even the $3$-spectral radii and $4$-spectral radii of $\beta$-stars with millions of vertices are gained with high probability and efficiency.

5.2 Approximation of Lagrangians of hypergraphs

When $p=1,$ the $1$-spectral radius is also known as the Lagrangian of a hypergraph (2.4). However, $f({\bf x})$ is not smooth at ${\bf x}$ who has some zero elements. We use $\lambda^{(p_{\vartheta})}(G)$ to approximate $\lambda^{(1)}(G)$, with $p_{\vartheta}$ being denoted as

[TABLE]

Since $\lim_{\vartheta\rightarrow\infty}p_{\vartheta}=1,$ we have $\lim_{\vartheta\rightarrow\infty}\lambda^{(p_{\vartheta})}(G)=\lambda^{(1)}(G)$ from Proposition 3.1. Therefore, we can use $p_{\vartheta}$-spectral radius to approximate the Lagrangian of a hypergraph. The function $f_{p_{\vartheta}}({\bf x})$ is continuous and differentiable and CSRH method is feasible for computing $p_{\vartheta}$-spectral radius of a uniform hypergraph. Let ${\bf w}$ be a vector such that its $i$th element being

[TABLE]

Then function $f_{p_{\vartheta}}({\bf x})=f_{p_{\vartheta}}({\bf w}^{[2\vartheta+1]})$ is also a semialgebraic function and satisfies the Łojasiewicz inequality (4.7). Therefore, the conclusions in Section 4 hold for $p_{\vartheta}$ in (5.3).

In this subsection, we show the results of CSRH method approximating Lagrangian of a hypergraph. First we give an example to demonstrate that the CSRH method is competent to compute the $p$-spectral radius of a uniform hypergraph when $p$ is a fraction in (5.3). Next, the numerical results of approximating the Lagrangians of complete hypergraphs by $p_{\vartheta}$-spectral radius are represented. The termination criteria of algorithms in the remaining part of this paper is set as $\|\nabla f({\bf x})\|\leq 10^{-6}.$

In Table 4, we present the consequences of the $p_{\vartheta}$-spectral radius of a 3-uniform $\beta$-star with 10 edges, with $p_{\vartheta}$ being the fraction in the first column. The true $p_{\vartheta}$-spectral radius can be acquired from Theorem 2.1. All experiments produce the exact $p_{\vartheta}$-spectral radius with probability $1$ and the relative error between our numerical result and the theoretical value obtained from Theorem 2.1 is at most $3.07\times 10^{-14}.$

An $r$-uniform hypergraph is said to be complete if it contains all possible edges when the number of its vertices is fixed. We use $C_{n}^{r}$ to denote a complete $r$-graph with $n$ vertices. Then the 3-graph $C_{4}^{3}$ is actually a tetrahedron with 6 edges. The Lagrangian of a complete uniform hypergraph can be obtained directly from Proposition 2.1.

We compute different $p_{\vartheta}$-spectral radii of 3 complete hypergraphs $C_{4}^{3}$, $C_{10}^{3}$ and $C_{20}^{3}$. In Figure 5, the ordinate reflects the error between the $p_{\vartheta}$-spectral radius and the true Lagrangian of the corresponding complete hypergraph which is obtained from the Proposition 2.1, while the abscissa means the value of $p_{\vartheta}-1.$ When $p_{\vartheta}$ approaches to $1,$ the $p_{\vartheta}$-spectral radius is close to the exact Lagrangian of the related hypergraph.

6 Network analysis

Not only the $p$-spectral radii, i.e., the optimal value of $f({\bf x})$ in (3.3), but also the optimal point ${\bf x}$ in (3.3) characterize the structure of hypergraphs. Recall (2.3) that an optimal point is called a $p$-optimal weighting. The elements of the $p$-optimal weighting reflect the importance of the corresponding vertices in the hypergraph. Therefore, we may call the $i$th element of the $p$-optimal weighting the impact factor of the $i$th vertex. Different selections of the parameter $p$ provide different criteria of the importance of the vertices. When $p$ is relatively large, the criterion tends to evaluate the importance of vertices more individually. When $p$ is relatively small, the ranking result demonstrates the significance of groups of vertices. In this section, we compute each $p$-spectral radius 10 times and choose the vector corresponding to the largest $f({\bf x})$ value as the $p$-optimal weighting.

6.1 A toy problem

We first employ a toy problem to illustrate the impact of the selections of $p$. We construct a 6-uniform weighted hypergraph with 8 edges as in Figure 6. The weights of all edges of this hypergraph are set as $1$, except the last one whose weight is $\frac{3}{2}$. Obviously from the hypergraph, the vertices numbered $1$, $31$, and $26$ are distinct from other vertices, and the edge $\{31,37,38,39,40,41\}$ is also distinct from other edges. In Table 5, we show the different ranking of vertices via different $p$-optimal weighting. The abbreviation Num. means the number of a vertex and Val. represents the impact factors of the corresponding vertices.

When $p=\frac{4}{3},$ the top $6$ vertices are in the edge who has the only largest weight among all edges. From Table 5, we can see that the impact factor of the top $6$ vertices in the $\frac{4}{3}$-optimal weighting are much greater than others. In fact, the value of all impact factors, except those corresponding to the top $6$ vertices, are less than $5\times 10^{-10},$ which means that the dominant vertices are the ones from the largest weighted edge and the others can be ignored. That is to say, the ranking in this case offers the most important group of the vertices. When $p=5$, the vertex numbered $26$ appears in the top 10 list and the difference among the top $10$ impact factors is not as great as that when $p=\frac{4}{3}$. When $p=16$, the top $3$ vertices are $1$, $31$, $26,$ and the impact factors of vertices that have same status in the hypergraph are rather close to each other. Then, we believe that the ranking results of $16$-spectral radius reflects the significance of vertices individually.

6.2 Author ranking

Ng et al. in [47] collected publication information from DBLP222http://www.informatik.uni-trier.de/ ley/db/ and gave different rankings of the authors according to different factors, such as citations of authors, category concepts, collaborations, and papers. In this subsection, we use the same data set in [47] and rank the authors based on their collaborations.333We would like to thank Dr. Xutao Li for providing the database.

We construct a weighted 3-uniform multi-hypergraph $G_{A}$ with $1,243,443$ edges to store the cooperation information. The vertex set is composed of numbers of the 10305 authors and each edge has 3 vertices indicating that these three authors have cooperations under a same topic. The weight of an edge is decided by the collaboration times among the three authors in this edge. The adjacency tensor of this multi-hypergraph $G_{A}$ is a sparse tensor with $1.17\%$ nonzero entries.

The example in Subsection $6.1$ shows that we can obtain the ranking score from different viewpoints by computing different $p$-optimal weighting. Therefore, we compute $2$-optimal weighting and $12$-optimal weighting of $G_{A}$ to get the author group ranking and the author ranking respectively. In Figure 7(a), the stars stand for the $2$-optimal impact factors of vertices of $G_{A}.$ Obviously, the majority elements of $2$-optimal weighting are extraordinarily close to zero and only dozens of corresponding stars are above the horizontal line of $y=0.1.$ In fact, $97.2\%$ of the entries in the $2$-optimal weighting are less than $10^{-3}$ and the elements that are greater than $0.1$ occupy only $1.8\%.$ On the other hand, the largest impact factor reaches to $0.4481$ and the upper stars are considerably larger than others. It means that the $2$-optimal weighting is dominated by a small proportion of its components and we regard these leading elements as a group. The top ten authors ranked according to the 2-optimal impact factor are presented in the second column of Table 6. The average collaboration times of each two authors among these top ten authors are $8.533,$ which is far larger than $9.76\times 10^{-4}$, the average collaboration times of each two authors among the whole $10305$ authors. Since these top ten authors have intimate cooperation, it is rational to consider them as a group and interpret the ranking in the second column as the most powerful group.

Stars in Figure 7(b) are the $12$-optimal impact factors of vertices of $G_{A}.$ The distribution of these stars is totally different from the ones in Figure 7(a). It can be seen in Figure 7(b) that the $12$-optimal impact factors of the $10305$ authors are uniform and most of them are concentrated in the internal between $0.006$ and $0.014.$ Because in the original data set, the collaboration times of different authors are mostly one or two and we rank the authors based on their collaborations, the balance and concentration of the impact factors match up with the cooperation information. The top ten authors generated via the $12$-optimal impact factors are listed in the third column of Table 6. Ng et al. also ranked the authors in the light of collaboration times and the influence of category concepts of their publications. We demonstrate the top 10 authors of their experimental result [47] in the MultiRank column in Table 6. It can be seen that $6$ of the top $10$ authors in the MultiRank are coincident with results of our $12$-optimal rank.

7 Conclusions

We convert the $p$-norm constraint in $p$-spectral radius problem into an orthogonal constraint, and propose a first order iterative algorithm CSRH for solving it. In this method, it is feasible to obtain a proper step length to satisfy the Wolfe conditions under the curvilinear line search. Convergence analysis shows that the CSRH method is globally convergent. The iterates converges to a $p$-optimal weighting. Numerical experiments show that CSRH method is efficient and powerful. In the author ranking application problem, we construct a weighted hypergraph with millions of edges. By computing $p$-spectral radius of this hypergraph, the most influential cooperation group and the top ten ranked authors are presented.