Implementing indexing on a SpectralRegion object (but first decide if we should!)

This is a follow-on from #790 (specifically the thread https://github.com/astropy/specutils/pull/790#discussion_r591854737). #790 implemented spectral indexing of the form ``spectrum[lowwl:highwl]``, using `SpectralRegion` as the underlying implementation. But what you *cannot* currently do is ``spectrum[SpectralRegion(...)]``. While it seems straightforward, @rosteen pointed out there's a subtlety around how `manipulation.extract_region` works at the edges vs how traditional python indexing works. So it's not clear what's the right way of interpreting the edges in the case of ``spectrum[SpectralRegion(...)]``.

This *may* mean this feature shouldn't get implemented at all, because we should just tell users to either use `extract_region` or "normal" indexing and not something sort of in-between. But it seems weird to me that the implementation is so close to ``spectrum[SpectralRegion(...)]`` but we can't do that directly:shrug:

Start with follow-on issue #790 and its linked discussion, then compare the existing spectrum[lowwl:highwl] behavior with manipulation.extract_region at region edges. The work is done when the edge semantics and whether SpectralRegion indexing should exist are explicitly decided, with the issue outcome reflected in the project behavior or guidance.

Written by the indexing model from the issue text.