Round 2 for ZSTD compression

[mobiusklein]

Describe the new term or terms you would like to add.

Revisits #167 and #168

Following on the discussion from the PSI Spring Workshop, since mzMLb did not receive uptake from the community for all the previously discussed reasons, we will revisit adding ZSTD compression to mzML.

The conversation proposes:

1. Provide ZSTD compression with byte shuffling with a "trade name".
2. Provide delta encoding followed by ZSTD compression with byte shuffling with a "trade name".
3. Do not propose separate ZSTD without byte shuffling because it introduces an extra complexity.

I ran some experiments, generate points from 100.0 to 2100.0 for double precision floats spaced by 0.01 using ZSTD compression level 9 and zlib compression level 9:

	Raw	Zlib	Shuffle + Zstd	Delta + Shuffle + Zstd
Bytes	1600000	570704	12040	560

I did not run a timing experiment because that will be language-sensitive.

Begin Bikeshedding

We wanted to use a trade name because it seemed to the room that a verbose name would be harder to understand.

Some randomly generated names:

1. ZSTD-Shuffle - The anti-trade name just call it what it is
2. MZSTD - The least creative adaptation of the domain
3. MZD - I was wrong, I could be less creative

[mobiusklein]

Some brief notes on complexity:

The byte shuffling algorithm is fairly straight-forward if you are willing to allocate a second array. In principle it is possible to implement a ping-pong shuffling algorithm that doesn't allocate but it will be very easy to corrupt memory with it.

The delta encoding algorithm is very easy to implement in-place or with allocation. I used the algorithm from mzMLb because the Numpress Linear algorithm does not really support byte shuffling (at least not the one I had access to).

[mobiusklein]

@chambm This conversation happened at the spring workshop, prompted by Sander Willems (@sander-willems-bruker) during the MzML for Ion Mobility and DIA session.

I'm willing to build this into ProteoWizard if you think this is acceptable. I believe there's already a static ZSTD library in the ProteWizard build tree so it shouldn't be a big change in maintenance burden.

[sander-willems-bruker]

Minor observation from my end. While I know zstd is an excellent compressor, it is also quite smart. The use case as presented here, i.e. perfectly even spaced points, is not a realistic one and esentially the "perfect" input for zstd. A test with random mz values in the range seems far more appropriate and will yield a less amazing but also more fair/realistic view of how zstd compared to zlib. Also note that I personally would recommend default levels for zstd. Yes you can tweak it to be more memory efficient but the extra time generally is not worth it.

[mobiusklein]

As is traditional, I took a random RAW file I had lying around on disk and ran some conversion benchmarks:

RB_161104_Janus_45_TiO2.raw from https://www.ebi.ac.uk/pride/archive/projects/PXD005931

5070 Profile MS1, 28026 Centroid MSn spectra

method	raw size	formatted	runtime (sec)
delta_zstd	1209986822	1.210 GB	17.4
zlib	1543244249	1.543 GB	100
zstd	1409583767	1.410 GB	17.3

In total, the MS1 m/z arrays from zstd used 85% of the space that they do in zlib and delta_zstd used 65% of the space of zlib. The ratio for MSn spectra was 88% and 77% respectively. Runtimes measured using the mzconvert example from https://github.com/mobiusklein/mzdata/tree/feature/zstd_compression with the thermo and zstd features enabled.

Compression ratios for MS1

[edeutsch]

This does looks quite encouraging, I hope it would not be difficult to implement in most applications.

[mobiusklein]

Following up. I implemented this in C++ as well, reminding myself why I have to use CMake the hard way. What other languages do we need to provide prototypes for? Java, Python and C#?

[chambm]

It's certainly easier to implement a new compression algorithm than a new HDF5-based format (even if all the metadata is mzML). But I don't know why we would expect much of a different trajectory for zstd (especially if it includes delta) as we saw for NumpressMS. Anyway Joshua if you want to make the PR to ProteoWizard we'll certainly incorporate it! It's all the other mzML readers out there we should worry about.

[mobiusklein]

@chambm Right now, there's an option to use delta encoding, or not, leaving it up to the user to choose. The thing that is always used is byte shuffling -> ZSTD compression. While delta encoding does improve compression on profile mode m/z arrays considerably, the bigger deal is the speed up encoding/decoding. Even so, 15% compression improvement over zlib-ng is nice, especially considering zlib-ng is considerably faster than the stock zlib library used by many tools.

Additionally, I'm writing implementations in languages I don't use so that when I have to publicize this I can point and say "just take the code off the shelf here".

"Single header" C++ implementation mzd.hpp probably will depend heavily on build system.

@lazear had initially requested this prior to mzMLb's failure to launch. Sage's popularity might help others use get on board with it?

@jcharkow I remember you had contributed to the discussion on this as well. Do you have any feedback on either name or implementation details?

[lazear]

+1 for byte shuffling and zstd. Happy to implement

[jcharkow]

Thanks for looping me into the discussion. Results look promising, especially because the file is both smaller and runtime is reduced as well!

Just a few requests:
It would be nice to also see benchmarks with msNumpress (followed by zlib compression) as well since this is something that I regularly use. Possibly worth also investigating msNumpress followed by zStd compression.

Also it would be nice to see some benchmarks on Bruker data since this is what @sander-willems-bruker had in mind when he originally proposed at HUPO-PSI (and that's what I primarily analyze). Any dataset you have access to is fine. If you need to download, I would suggest a single cell dataset as they tend to be smaller. (e.g. https://massive.ucsd.edu/ProteoSAFe/dataset_files.jsp?task=7c2180ca9eff4640b60576aa500d8721#%7B%22table_sort_history%22%3A%22main.collection_asc%22%7D)

On another note, Hannes also suggested investigating mspack, however, I am not sure how relevant this is because it does not support random access.

[mobiusklein]

@jcharkow I have read the mspack paper (10.1093/bioinformatics/btab636). This appears to combine out-of-band data that is compressed and standard gzip compression of the XML metadata. This is pretty much the same idea mzMLb, just that HDF5's chunked compression would allow random access on the XML metadata too.

I fetched E_20221108_EvoSepOne_3rdGenAurora15cm_CC_40SPD_whisper_scLF1108_M10_S1-A1_1_3103.d from ftp://massive-ftp.ucsd.edu/v06/MSV000093867/.

The whole .d directory is 4.06GB, but the analysis.tdf_bin file is only 1.92GB, and a second file linespectra is also 1.94GB though I do not know its purpose. Let's call the core of the data a nice even 2GB since the analysis.tdf is 82MB.

I modified the mzconvert code so it wouldn't do some custom signal processing and ran it through

method	raw size	formatted	runtime (sec)
zlib	6552102347	6.6GB	1551.37
zstd	6875499979	6.9GB	472.50
delta-zstd	5928218176	5.9GB	488.54

A note about compression levels, I use the "best" compression level for zlib-ng, while using the "default" compression level for zstd.

I will update when I finish testing the MS-Numpress method. Because https://crates.io/crates/numpress only supports doing linear encoding and only on 64-bit floats, I can't do a complete comparison with slof compression for the intensities.

[jcharkow]

@jcharkow I have read the mspack paper (10.1093/bioinformatics/btab636). This appears to combine out-of-band data that is compressed and standard gzip compression of the XML metadata. This is pretty much the same idea mzMLb, just that HDF5's chunked compression would allow random access on the XML metadata too.

I fetched E_20221108_EvoSepOne_3rdGenAurora15cm_CC_40SPD_whisper_scLF1108_M10_S1-A1_1_3103.d from ftp://massive-ftp.ucsd.edu/v06/MSV000093867/.

The whole .d directory is 4.06GB, but the analysis.tdf_bin file is only 1.92GB, and a second file linespectra is also 1.94GB though I do not know its purpose. Let's call the core of the data a nice even 2GB since the analysis.tdf is 82MB.

I modified the mzconvert code so it wouldn't do some custom signal processing and ran it through

method raw size formatted runtime (sec)
zlib 6552102347 6.6GB 1551.37
zstd 6875499979 6.9GB 472.50
delta-zstd 5928218176 5.9GB 488.54
A note about compression levels, I use the "best" compression level for zlib-ng, while using the "default" compression level for zstd.

I will update when I finish testing the MS-Numpress method. Because https://crates.io/crates/numpress only supports doing linear encoding and only on 64-bit floats, I can't do a complete comparison with slof compression for the intensities.

Thank you for the comments and for doing these benchmarks! I will also convert this file with diaPysef (which uses msnumpress + zlib compression) for comparison. Runtime will be less relevant here since it is python so it is much slower.

[jcharkow]

I've done a benchmark with msnumpress on my end for the same file that you benchmarked with above.

m/z array is 64-bit float with MS-Numpress linear prediction compression followed by zlib compression
intensity array is 64-bit float with MS-Numpress short logged float compression followed by zlib compression
ion mobility array is MS-Numpress short logged float compression followed by zlib compression.

The file size is an even 5.0Gb which is smaller than the files above. However msNumpress is lossy so that might not be ok in some contexts, even though it works for me well here.

Conversion takes 25 minutes but a lot of this is likely to do with the fact that was done in python (using diaPysef )

@mobiusklein do you mind also running delta-zstd on a higher compression level? I am curious with how this compares to the filesize of my MSNumpress results.

Also I wanted to double check that you were also compressing the ion mobility array.

Thanks!

[mobiusklein]

It has been a while since I looked in detail at Numpress, I forgot that it stores data not just in an easier-to-compress format, but also using a lower width data type, with 64 bit floats becoming 32 bit integers. This makes the comparison with the delta encoding uneven since it encodes the deltas as 64 bit floats, and then asks the compressor to collapse them.

method	raw size	formatted	runtime (sec)
zlib	6552102347	6.6GB	1551.37
zstd	6875499979	6.9GB	472.50
delta-zstd	5928218176	5.9GB	488.54
numpress	5093744855	5.1GB	1219.59
numpress-zlib	4838844933	4.8GB	1312.70

This used the MS-Numpress linear encoding on the m/z arrays, and zlib on the ion mobility array. It turned out that in all of these benchmarks used zlib on the ion mobility arrays.

Worse, I found that when you apply the byte shuffling + zstd compression, ion mobility compression is much, much worse than zlib or zstd on its own, and zstd compression isn't competitive using the default compression level. Even more problematic is that not every spectrum's ion mobility array achieves better compression levels with higher zstd compression levels. What is surprising is that zlib is able to compete

To some extent, this is reasonable. The ion mobility array is not sorted, so there is less structure to derive patterns from. We also don't see consistent improvements in compression for zstd on intensity arrays.

Here's a plot of compression ratios relative to the default level per spectrum for this file:

Values below 1.0 are spectra where the lower compression level actually performs better than the higher level.

This, combined with the fact that ion mobility instrumentation is more coarsely digitized compared to m/z, suggests that we should be using a less "clever" algorithm for compressing profile ion mobility data. We know we have a near infinite # of possible m/z values, but a fixed # of ion mobility values which we can identify with a 2 byte integer. In principle, we could do very simple dictionary encoding. I will follow up when I finish doing a little bit of experimentation.