Can PC1 explain more than 90% of variance? I've been running some analysis using morphological data. However, when I run a PCA, PC1 explains more than 90% of variance. I'm not sure if this is possible. If not, what should I need to do?

 A: In case you need a non-simulated example: I did help.search("morpholog", agrep = FALSE) to find morphological data sets in the set of R packages I happen to have installed on my system.  A useful one is ade4::tortues, a data set on 48 painted turtles from Jolicoeur and Mosimann 1960.
long, larg, and haut are the length, width, and height of the turtles in mm; because they're all measured in the same units, it shouldn't be necessary to scale the variable before/while computing principal components.
library(ade4)
## column 4 is sex
pairs(tortues[,1:3], col = tortues[,4], gap = 0)
summary(prcomp(tortues[,1:3]))

PC1 has 98% of the variance. (Jolicoeur and Mosimann get 97.61%; I get 97.16% if I set scale= TRUE; haven't checked for other issues [transcription errors etc.]).
Importance of components:
                           PC1     PC2     PC3
Standard deviation     25.3100 2.40272 2.26449
Proportion of Variance  0.9833 0.00886 0.00787
Cumulative Proportion   0.9833 0.99213 1.00000


Jolicoeur, P. and Mosimann, J. E. (1960) Size and shape variation
in the painted turtle. A principal component analysis. Growth,
24, 339-354.
A: Let $S \in \mathbb{R}^{p \times p}$ be the sample covariance matrix (which equals to $X'X$ for centered and intercept-free design matrix $X$), and $\lambda_1 \geq \lambda_2 \geq \cdots \geq \lambda_p \geq 0$ be all the non-negative eigenvalues of $S$, then the total variation percentage explained by the first $k$ PCs is given by
\begin{align}
p_k = \frac{\sum_{j = 1}^k \lambda_j}{\sum_{j = 1}^p \lambda_j}, \; k = 1, 2, \ldots, p.
\end{align}
In particular, the total variation percentage explained by the first PC is
\begin{align}
p_1 = \frac{\lambda_1}{\sum_{j = 1}^p\lambda_j}. \tag{1}
\end{align}
In view of $(1)$, theoretically there are no upper bounds that are strictly less than $1$ for $p_1$.  In other words, there exists $X$ such that $p_1$ can be arbitrarily close to $1$. For example, any design matrix $X$ in the form of $(2)$ below satisfies that its first PC explains $1 - \varepsilon$ total variation:
\begin{align}
X = U\begin{bmatrix}
\operatorname{diag}(\sqrt{1 - \varepsilon}, c, \ldots, c) \\
0 
\end{bmatrix}V', \tag{2}
\end{align}
where $c = \frac{\sqrt{\varepsilon}}{\sqrt{p - 1}}$, $U$ and $V$ are arbitrary order $n$ and order $p$ orthogonal matrices. Therefore, it is not surprising at all to observe cases such that $p_1 > 90\%$.
A: Certainly that can happen. Below is a simulated example in R with just two original variables. As long as they are strongly enough correlated, they can very well be summarized by the projection to a single line, which means that the first principal component explains an arbitrarily large proportion of total variance:
library(MASS)
set.seed(1)
dataset <- mvrnorm(1000,c(0,0),cbind(c(1,0.99),c(0.99,1)))
summary(prcomp(dataset))
plot(dataset,las=1,pch=19,cex=0.6)

Output:
Importance of components:
                         PC1     PC2
Standard deviation     1.460 0.10400
Proportion of Variance 0.995 0.00505
Cumulative Proportion  0.995 1.00000


I very much recommend Making sense of principal component analysis, eigenvectors & eigenvalues.
