How do you know what the correct Bode plot is that your taking as your reference for comparison?
Basically repeating from Star's comment for completeness: In Matlab, you are geting the Bode plot from the second input to the first output. Is that what you want? From the name s1a1 and from the title on the Mathematica plot, it looks like you may have wanted the Bode plot from the first input to first output. And be careful with the Hz vs. rad/s thing.
The bode function in the CST takes great care in selecting its frequency vector if not specified by the user (as in your case) in an attempt to make sure it captures important dyanmics. In your case, all of the poles and zeros of s1a1 are basically on the imaginary axis. If you don't catch exactly the right frequency to evaluate you'll miss a peak. Try this and see what you get:
What frequency vector did Mathematica use?
Finally, once you have the state space model, there's no reason convert to tf first. In fact, I'm quite certain that documentation for older versions of the CST specifically recommended NOT doing this. Don't know if that's the case now, but I'm pretty sure that if you start with a ss representation, you should just use it and not bother converting to the other forms. If all you wanted is the Bode plot from the second intput to the first output you get that directly: