This is caller "curvature" and has been one of the main
challenges since the first bipolar bandgap designs. Many
papers on curvature correction out there, and many of
them are blind alleys you don't want to go down (for
example, the idea that you should "correct" the 2nd
order error with a cubic term tends to bring more pain
than relief).
The parabolic character sets a minimum deflection
over temp which limits your ability to trim to a low
across-envelope max error.
Base current and Rb is one error source. So is the
higher order tempco (2nd, 3rd, ...) of resistors. A
dead flat metallic resistor does not suffice to eliminate
it all, but at least gets rid of the resistor contribution
(if you can have it - used to be a common feature
of older precision analog processes, but seldom seen
in digital or "mixed signal" (i.e. digital plus MIM cap)
foundry flows.
I suspect plenty of process engineering at ADI and
Linear Tech for the resistors of the special flows their
high precision continuous-time reference products run
on. We were never able to come close to their art
using more generic linear processes, even ones with
thin film metallic resistors.
Another approach involves a digital cal-map, a temp
sensor, a crude ADC and a trim-DAC pushing a
compensating current into the bandgap core, or
a voltage summed with its output. In higher density
nodes this can be more area-economical than an
E-trim (fuse, laser) block. It can take you to a dead
flat tempco at the output.