How do I find the return loss of a transformer from its 4-port S-parameter data?
By reading some application notes, I discovered, indeed, that the 4-port S-parameters measurements file is meant for creating a model of the transformer, and it can be used with a simulator directly, such as Genesys or ADS. Those are expensive tools unavailable for hobbyists, but I also found the free and open source Qucs simulation supports n-port S-parameter file as well.



To use the s4p file for simulation...


Create two power sources and specify their input impedances, in this case, 50 Ω and 800 Ω, and ground them.
Create two equations for S11_dB and S21_dB, and create a "S-parameter simulation", set the simulation type to log, set the start to stop frequency, and set the number of points.


These steps are needed for all simulations. A quick way of setting the skeleton up is using the filter synthesis wizard to obtain a random filter, and delete all the R, L, C components. 


Add a "n-port S-parameter file". edit the properties and change the number of port to "4", set the file location to s4p file provided by Mini-Circuit. Since port 1/2 are the primary side, port 3/4 are the secondary, we need to rotate the device by 90 degrees, a bit awkward but it does the job, finally ground the Ref port. 


Run the simulation.



It works, we get identical results just like the datasheet (the only difference is the definition of positive/negative)! We can now even plot the results in different format, such as a Smith chart and a VSWR chart, hooray!

There is still ~4 dB of difference, probably due to measurements under different conditions, or due to artifacts of the simulator or its setup, that worth further investigation. 

Update: Agilent ADS produces the same simulation result. Thus it's probably related to measurement conditions.