magnetism is a result of the virtual photons interactions because electromagnetism is a result of it.
Now why we have electromagnetism, again I am going to address you to a special relativity book, or Griffiths introduction to electrodynamics - because there he is trying to show mainly with words how magnetism is a result of the electrism depending on the reference frame of the observer. The mathematical result of his reasoning is the ability/existence of a 2rank antisymmetric tensor which contains the Electric and Magnetic field.
Another way to see how those things are connected, is via the Maxwell's equations- there you can see how the electric field can really generate magnetic field and vice versa. Of course everything becomes much simpler and clear in special relativity.
In fact there's a joking comment from most people who teach the Maxwell equations, that what Special Relativity intended to do, was to make the Maxwell Equations frame invariant (they are not Gallilean invariant), as can be easily checked. I hope this thread helps somewhat:
https://www.physicsforums.com/showthread.php?t=615937
The main problem with Gallilean transformations is that they don't change the time...
Nevertheless... Once you are able to work in special relativity, you are able as a said to define the antisymmetric 2rank tensor (field strength tensor field):
F_{\mu \nu} = - F_{\nu \mu}
F= dA
or in components:
F_{\mu \nu} = \partial_{[\mu} A_{\nu]} \equiv \partial_{\mu} A_{\nu} -\partial_{\nu} A_{\mu}
Where A^{\mu}= (\Phi, \vec{A})
is the electromagnetic 4 potential. The zero-th component is just the electric potential field, where the spatial components is the magnetic potential vector field.
There are many ways to define the above F_{\mu \nu}. One way, in particular in particle physics, is via the Invariance of the lagrangian under local U(1) transformations, which accept such terms. This is a nice way because it's also able to be extended for other groups (non-abelian... U(1) is an abelian group) and thus even give you the corresponding field strength tensors for SU(3) (QCD), or SU(2) (somewhat the Weak Interactions- of course not exactly since you need to learn Higgs Mechanism for it). Nevertheless, I shall continue with the particle way of seeing things. The Lagrangian for a fermion for example, is just the Dirac Lagrangian:
L= i \bar{\Psi} \gamma^\mu \partial_{\mu} \Psi
If you try to use a local U(1) transformation, to the Dirac spinors, which is equivelant to making the change:
\Psi \rightarrow e^{ig a(x)} \Psi
(a global U(1) is in case a doesn't depend on spacetime points x^{\mu}...In that case your lagrangian is invariant as it is).
You can see that the above Lagrangian is not invariant... it will change. In order for it not to change you have to allow a minimal coupling- that means to change the way you do partial derivatives... The same thing you do with General Relativity, where instead of the normal derivatives you have to allow a minimal coupling to the gravitational field through the Christoffel symbols/connections.
It happens to be a correct way to allow:
\partial_{\mu} \rightarrow \partial_{\mu} + ig A_{\mu} \equiv D_{\mu}
After replacing the partial derivatives with D_{\mu} the Lagrangian is able to become invariant, if you allow A_{\mu} to accept a U(1) transformation:
A_{\mu} \rightarrow A_{\mu}'= A_{\mu} + G \partial_{\mu} a(x)
G some factor...
The last is also an indication that your field is somewhat a photon - it accepts the same transformations you can find for a photon in special relativity. By the above, also a local U(1) invariant term in the Lagrangian is the quantity: F^2 = F_{\mu \nu} F^{\mu \nu}
This gives the photon's kinetic terms..
L_{photon kin}= - \frac{1}{4} F^2
Why? because by making a change of the A_{\mu} as above, and because of the antisymmetry you will have:
F_{\mu \nu}' = F_{\mu \nu} + G \partial_{\mu} \partial_{\nu}a(x)-G \partial_{\nu} \partial_{\mu}a(x)= F_{\mu \nu}
All the above might seem like magic, that we are able to reconstruct the electromagnetic tensor from imposing a local U(1) transformation. Also most of the points can be exctracted only by looking at the topological properties of U(1) gauge group (seems like geometry), something you can find in Peskin & Schroeder Chapter 15 ( particularly paragraph 1). That's why some people mentioned that Electromagnetism is in fact a result of gauge invariance.Nevertheless... The next step for someone is to see how fields interact with the photons... In fact this is more difficult for me to explain, because the vollume of what must be said is huge and can result in writing a whole book- which I am unable to do at the moment... I already wrote too much in order to help you see how this is connected to gauge symmetry group.
In order for someone though to check how it results, is to write down the full lagrangian - allowing all the fields that you want, and allowing them to interact somehow. The charged fields will act as a source of electromagnetic field (you can check that by looking at the equations of motion for F_{\mu \nu} (in fact in the easiest way when you have no interactions of photons with fields, the equations of motion for F are going to give you the Maxwell equations. Some clever mind used this fact in order to create a t-shirt with:
"And god said F=dA and so there was light"
much more simple than writing the whole 4-set of Maxwell's equations, isn't it?

In F's equations of motion alone, are not enough to allow both the 4 equations to be reconstructed. To get them you need either to allow the dual F_{Dual}^{\mu \nu} = \epsilon^{\mu \nu \rho \sigma} F_{\rho \sigma} or dF=0 \rightarrow \partial_{[\mu} F_{\nu \rho]}=0 which you can see it holds because dF= d(dA)=0.
\epsilon is the Levi Civita symbol in 4 dimensions (or totally antisymmetric element).
\epsilon^{\mu \nu \rho \sigma}= +1,~~ (\mu \nu \rho \sigma) even permutations of (0123) cyclic group.
\epsilon^{\mu \nu \rho \sigma}= -1,~~ (\mu \nu \rho \sigma) odd permutations of (0123) cyclic group
\epsilon^{\mu \nu \rho \sigma}= 0,~~ else (in fact if any index is the same as other)
eg
\epsilon^{0123}=+1, ~~ \epsilon^{1032}=+1, ~~ \epsilon^{1023}=-1, ~~\epsilon^{1302}=-1, ~~ \epsilon^{0120}=0,~~ \epsilon^{1013}=0
etc
Nevertheless, once you have the full lagrangian, you can start looking at Feynman Diagrams. The Lagrangian will indicate you what couplings you can have in a Feynman Diagram, and from that you can start looking at charged particle interactions with the Electromagnetic field - there you will need the photons as propagators (virtual/off-shell particles) and thus you see how particles interact with those virtual photons.
That's the electromagnetic interactions- in them you have hidden both the electric and magnetic interactions, but they are identical : two views of the same coin.
I hope this is enough to motivate you into getting to study things like differential geometry, group theory and quantum field theory.