We present the results on topological susceptibility and chiral observables in $N_f=2+1+1$ QCD for temperature range $120<T<600$ MeV. The lattice simulations are performed with Wilson twisted mass fermions at physical pion, strange and charm masses.
In high-$T$ region $T\gtrsim 300$ MeV the chiral observables are shown to follow leading order Griffith analyticity, and the topological susceptibility follows a power-law decay as in the instanton dilute gas models. The measured topological susceptibility is used to estimate the mass of QCD axion. The resulting axion mass constraints are in agreement with our previous studies at higher pion masses.