Dear Mattia,
I am unsure if there's a more direct way but if that has not been implemented, the following way, which is not smart, might help.
one way to do it is:
1. for each site, create a model in the format of nexus where the exchangeability matrix is the one you are going to use, say LG, and the equl freq vector is the site-specific freq obtained from .sitefreq [put it at the last line in the nexus file (see PAML's "dat/.dat")]
2. use "-mdef model.nex -m model" to import the model
3. iterate the above for each site
one thing you might want to be careful about is that by giving a new Freq the exchangeability matrix R=(r_ij) need to be normalized but i think iqtree will automatically adjust it for you so hopefully no need to worry.
an alternative way is to use
iqtree -m LG+F{F{0.123119/0.019475/0.019852/0.018583/0.018711/0.017275/0.018723/0.050388/0.016402/0.119697/0.161399/0.012776/0.035838/0.057019/0.030913/0.043472/0.049935/0.012600/0.039929/0.133894}
(remember to replace the freq w/ your site freq and do the same for each site)
kindly note that I have tested neither of the above but my gut feeling is that they both work, although not in a very easy way. thx!
best,
sishuo