Summary
Below, we prove that
$$
\|\hat\Phi\|_{\rm op}=\sup_{\rho\in D(\mathcal{X})}\sqrt\frac{\gamma(\Phi(\rho))}{\gamma(\rho)}
$$
where $D(\mathcal{X})$ denotes the set of density matrices on the Hilbert space $\mathcal{X}$ associated with the input of the quantum channel $\Phi$ and $\gamma(\rho)=\mathrm{tr}(\rho^2)$ is the purity of $\rho$.
This justifies the interpretation of $\|\Phi\|_{\rm op}$ as the square root of the largest multiplicative increase in purity due to the action of the channel $\Phi$. As a corollary, we obtain $\|\Phi\|_{\rm op}\in[1,\sqrt{d}]$ where $d$ is the dimension of the input Hilbert space $\mathcal{X}$.
Channel operator norm
Operator norm $\|.\|_{\rm op}$ is defined$^1$ as
$$
\|A\|_{\rm op}=\sup_{\|v\|_2\le 1}\|Av\|_2=\sup_{v\ne 0}\frac{\|Av\|_2}{\|v\|_2}.\tag1
$$
where $\|v\|^2_2=\sum_k |v_k|^2$ is the vector $\ell_2$ norm on the domain and codomain of $A$.
If $\hat\Phi$ is the matrix representing the quantum channel $\Phi:L(\mathcal{X})\to L(\mathcal{Y})$, then for all $X\in L(\mathcal{X})$
$$
\mathrm{vec}({\Phi(X)}) = \hat\Phi\,\mathrm{vec}(X)\tag2
$$
where $\mathrm{vec}$ is the vectorization map. Moreover, the $\ell_2$ norm of the vectorization of $X$ coincides with the Frobenius norm $\|.\|_F$ of $X$
$$
\|\mathrm{vec}(X)\|_2=\sqrt{\sum_{i,j=1}^n|x_{ij}|^2}=\|X\|_F\tag3
$$
where $n=\dim\mathcal{X}$. Thus, from $(1)$ we have
$$
\|\hat\Phi\|_{\rm op}=\sup_{\|X\|_F\le 1}\|\Phi(X)\|_F=\sup_{X\ne 0}\frac{\|\Phi(X)\|_F}{\|X\|_F}.\tag4
$$
Matrix achieving supremum is positive semidefinite
We will prove that the supremum in $(4)$ is achieved for a positive semidefinite $X$. First note that $\|\Phi(.)\|_F$ is continuous and $\{X\,|\,\|X\|_F\le 1\}$ is compact, so the supremum is achieved for some operator $M$. By linearity of $\Phi$ and absolute homogeneity$^2$ of the norm, $\|M\|_F=1$. We can write $M$ as the sum of its Hermitian and anti-Hermitian parts
$$
M=\frac{M+M^\dagger}{2}+\frac{M-M^\dagger}{2}.\tag5
$$
Setting $\alpha:=\frac12\|M+M^\dagger\|_F$ and $\beta:=\frac12\|M-M^\dagger\|_F$ we can rewrite $(5)$ as
$$
M=\alpha H+i\beta G\tag6
$$
where $H=\frac{M+M^\dagger}{2\alpha}$ and $G=\frac{M-M^\dagger}{2i\beta}$ are Hermitian matrices with unit Frobenius norm. Moreover, for any Hermitian $A,B$ we have $\|A+iB\|_F^2=\|A\|_F^2+\|B\|_F^2$, so $\alpha^2+\beta^2=1$.
Now, $\Phi$ sends Hermitian matrices to Hermitian matrices, so
$$
\|\Phi(M)\|_F^2=\alpha^2\|\Phi(H)\|_F^2+\beta^2\|\Phi(G)\|_F^2.\tag7
$$
But $\|\Phi(.)\|_F^2$ is convex$^3$, so $\|\Phi(M)\|_F^2\le\max\{\|\Phi(H)\|_F^2, \|\Phi(G)\|_F^2\}$. Therefore, we can take $M$ to be Hermitian.
Writing $M=R-S$ for some positive semidefinite $R$ and $S$ with orthogonal supports, we have $|M|=R+S$ and
$$
\begin{align}
\|\Phi(M)\|_F^2&=\|\Phi(R)-\Phi(S)\|_F^2\\
&=\|\Phi(R)\|_F^2+\|\Phi(S)\|_F^2-2\langle\Phi(R),\Phi(S)\rangle_{HS}\\
&\le\|\Phi(R)\|_F^2+\|\Phi(S)\|_F^2+2\langle\Phi(R),\Phi(S)\rangle_{HS}\\
&=\|\Phi(R)+\Phi(S)\|_F^2\\
&=\|\Phi(|M|)\|_F^2
\end{align}\tag8
$$
where we used the fact that the Hilbert-Schmidt inner product of two positive semidefinite operators is a non-negative real number. Therefore, the supremum in $(4)$ is achieved for a positive semidefinite $X$ and we can rewrite $(4)$ as
$$
\|\hat\Phi\|_{\rm op}=\sup_{\rho\in D(\mathcal{X})}\frac{\|\Phi(\rho)\|_F}{\|\rho\|_F}\tag9
$$
where $D(\mathcal{X})$ denotes the set of density matrices on Hilbert space $\mathcal{X}$.
Relation to purity
Finally, Frobenius norm $\|\rho\|_F$ of a density matrix $\rho$ is the square root of its purity $\gamma(\rho)$. Therefore, $(9)$ becomes
$$
\|\hat\Phi\|_{\rm op}=\sup_{\rho\in D(\mathcal{X})}\sqrt\frac{\gamma(\Phi(\rho))}{\gamma(\rho)}.\tag{10}
$$
This confirms the interpretation of $\|\hat\Phi\|_{\rm op}$ as the square root of the maximum multiplicative increase in purity due to the action of the channel $\Phi$.
In $d$-dimensional Hilbert space we have $\gamma(\rho)\in [\frac{1}{d},1]$, so $\|\Phi\|_{\rm op}\le\sqrt{d}$. Also, by Schauder's theorem, $\Phi$ has a fixed point, so $\|\Phi\|_{\rm op}\in[1,\sqrt{d}]$.
$^1$ Strictly speaking, operator norm may be induced by any pair of vector norms on the domain and codomain. Since operator norm $\|A\|_{\rm op}$ induced by the $\ell_2$ norm $\|v\|^2_2=\sum_k |v_k|^2$ is equal to $A$'s largest singular value, we assume that this is the operator norm the question is concerned with.
$^2$ Absolute homogeneity means that $\|sX\|_F=|s|\|X\|_F$ for all $s\in\mathbb{C}$. As pointed out by @Markus Heinrich in the comments, we may alternatively appeal to the convexity of $\|\Phi(.)\|_F$ and $\{X\,|\,\|X\|_F\le 1\}$.
$^3$ Because $\Phi$ is linear, $\|.\|_F$ is convex and $f(x)=x^2$ is convex and non-decreasing on $[0,+\infty)$.
