Depending on what the data look like, it may be much easiers, quicker, and more efficient to just use sort() or maxk() rather than findpeaks().
Power spectral density (computed by pwelch()) typically results in large spike rather than broad curves so the two largest spikes would merely me the two largest data points.
[pks,lcs] = sort(y_spectrum,'descend');
The two tallest peaks are located at indices lcs(1:2) and their peak values are pks(1:2).
Alternatively,
[pks,lcs] = maxk(y_spectrum;
If the data are not spikes then Star Strider's answer below is indeed safer. Here's a demo comparing the two methods. Note that maxk() method is 84x faster than findpeaks()+maxk(). maxk() method is also 15x faster than sort().
% Load built-in data
load handel.mat
y=y';
y_spectrum=pwelch(y);
% maxk method
n = 2; %find n tallest peaks
[pks,lcs] = maxk(y_spectrum,n);
%findpeaks + maxk method
[pks2, lcs2]=findpeaks(y_spectrum);
[maxpks,idx] = maxk(pks2, n);
max2lcs = lcs2(idx);
% Plot results
figure();
plot(y_spectrum,'k-','DisplayName','PSD')
hold on
plot(lcs,pks,'r*','DisplayName', 'maxk')
plot(lcs2(idx), maxpks, 'bs','DisplayName','findpeaks+maxk')
legend()
