Markov networks are frequently used in sciences to represent conditional independence relationships underlying observed variables arising from a complex system. It is often of interest to understand how an underlying network differs between two conditions. In this paper, we develop methods for comparing a pair of high-dimensional Markov networks where we allow the number of observed variables to increase with the sample sizes. By taking the density ratio approach, we are able to learn the network difference directly and avoid estimating the individual graphs. Our methods are thus applicable even when the individual networks are dense as long as their difference is sparse. We prove finite-sample Gaussian approximation error bounds for the estimator we construct under significantly weaker assumptions than are typically required for model selection consistency. Furthermore, we propose bootstrap procedures for estimating quantiles of a max-type statistics based on our estimator, and show how they can be used to test the equality of two Markov networks or construct simultaneous confidence intervals. The performance of our methods is demonstrated through extensive simulations. The scientific usefulness is illustrated with an analysis of a new fMRI data set.