title: Long term optical and infrared variability characteristics of Fermi Blazars
authors: P. Z. Safna C. S. Stalin Suvendu Rakshit Blesson Mathew
text: 
  
  
  firstpage–lastpage
  Statistical properties of eigenvalues of an ensemble of pseudo-Hermitian Gaussian matrices
      M. P. Pato
      December 30, 2023
  ==========================================================================================
  
  
  
  
  
  We present long term optical and near infrared flux variability analysis of 37 blazars detected in the γ-ray band by the Fermi Gamma-Ray Space Telescope. Among them, 30 are flat spectrum radio quasars (FSRQs) and 7 are BL Lac objects (BL Lacs). The photometric data in the optical (BVR) and infrared (JK) bands were from the Small and Moderate Aperture Research Telescope System acquired between 2008-2018. From cross-correlation analysis of the  light curves at different wavelengths, we did not find significant time delays between variations at different wavelengths, except for three sources, namely PKS 1144-379, PKS B1424-418 and 3C 273. For the blazars with both B and J-band data, we found that  in a majority of FSRQs and BL Lacs, the amplitude of variability (σ_m) in the J-band is larger than B-band consistent with the dominance of the non-thermal jet over the thermal accretion disc component. Considering FSRQs and BL Lacs as a sample, there are indications of σ_m to increase gradually towards longer wavelengths in both, however, found to be statistically significant only between
  B and J-bands in FSRQs. In the B-J v/s J colour magnitude diagram, we noticed complicated spectral variability patterns. Most of the objects showed a redder when brighter (RWB) behaviour. Few objects showed a bluer when brighter (BWB) trend, while in some objects both BWB and RWB behaviours were noticed. These results on flux and colour characteristics indicate that the jet emission of FSRQs and BL Lacs is indistinguishable.
  
  
  
  
  
  galaxies:active – (galaxies:) BL Lacertae objects: general – galaxies:nuclei – galaxies:jets
  
  
  
  
  
  
  
  § INTRODUCTION
  
  
  Active Galactic Nuclei (AGN) are actively accreting super massive black holes (SMBH) at the centers of massive galaxies <cit.>. The radiation from AGN is believed to result from the accretion of matter by the SMBH at the center of its host galaxy <cit.>. Among AGN are the blazars that show 
  large amplitude flux variations <cit.> over the entire accessible bands of the electromagnetic spectrum. The emission from blazars is dominated by non-thermal emission from their relativistic jets that are aligned close to the observer. They are broadly classified into two categories namely flat spectrum radio quasars (FSRQs) and BL Lac objects (BL Lacs). This division
  is from an observational perspective and it is based on the presence of broad emission lines in their optical
  spectra with FSRQs having broad emission lines with equivalent width (EW) > 5 Å. On a more physical ground, it is now thought that FSRQs have radiatively efficient optically thick and geometrically thin accretion disk, while BL Lacs have radiatively inefficient accretion disk that cannot photoionize the 
  line emitting clouds <cit.>. 
  
  The broad band spectral energy distribution (SED) of blazars
  show a two hump structure. The low energy hump, extending from radio through X-rays is dominated by the optically thin synchrotron emission process by the energetic 
  electrons in the jet <cit.>. The radio/optical emission that contributes to this component of the SED is found to be highly polarized <cit.>. The origin of the high energy hump that spans the X-ray and γ-ray region is not well understood. In the leptonic scenario <cit.>, the high energy hump could be produced  by inverse Compton scattering of seed photons by the same electrons 
  in the jet that are responsible for the synchrotron emission (synchrotron self Compton; ) as well as photons that are exterior to the jet (external Compton;).  However, the nature of the seed photons that take part in the inverse
  Compton process is poorly understood.  Alternatively, in the hadronic scenario, the second hump could also be produced by synchrotron radiation from protons that are co-accelerated along with the jet electrons, interaction of relativistic protons with
  external radiation fields as well as  proton induced particle 
  cascades <cit.>.
  Clues to a better understanding of the emission processes in blazars can be attained by studying correlated flux variations 
  between different spectral bands. 
  
  
  
  
  Blazars have been extensively studied for flux variability in multiple wavebands
  <cit.>, 
  yet, we do not have a clear understanding of the cause of flux variations in them. In addition to flux variations, blazars also exhibit spectral variations. Available observational results in literature point to complex spectral variability behaviours in 
  them <cit.>. Analysis of colour magnitude trends in monitoring data of blazars has revealed complex patterns such as bluer-when-brighter (BWB) trend <cit.>, 
  redder when brighter (RWB) trend <cit.>, both BWB and RWB trends <cit.> and no/weak spectral change with brightness <cit.>. To 
  better understand the complex flux and spectral behaviour of blazars it is important to analyse near simultaneous data on a large number of blazars that spans over many years. With this objective we have carried out a systematic analysis of the photometric monitoring data on a total of 37 blazars to characterize their flux and spectral variability behaviour as well as correlations if any between flux variations between different wavelengths. The sample  and data are described in Section 2. The analysis is discussed in Section 3. The results are discussed in Section 4 followed by the summary in the final section. 
  
  
  
  § SAMPLE AND DATA 
  
  Our sample of sources for this study was selected from the list of blazars that 
  are routinely monitored by the Small and Moderate Aperture Research 
  Telescope System (SMARTS[http://www.astro.yale.edu/smarts/fermi], )
  in the optical and near infrared bands. They are part of a list of γ-ray 
  emitting blazars regularly monitored by SMARTS. Firstly, we selected all objects from SMARTS 
  and restricted our analysis to only those objects that satisfy the criteria of 
  the availability of greater than 40 data points in at least any two optical 
  bands and one infrared band during the period 2008 - 2018.  
  This constraint leads us to a total  of 37 blazars of which
  30 are FSRQs and 7 are BL Lacs. For 19 sources 
  in our sample, we have data in all the five photometric bands namely B, V, R, J and K. For 2 sources we have
  data in B, V, R and J bands. Data in three bands B, R and J are available for 12 sources in our
  sample and the remaining 4 sources have data in four bands, namely B, R, J and K.
  The details of the objects that were analysed for flux variations in this work are given in Table <ref>. 
  From visual inspection of each of the light curves, we noticed few photometric points with large error bars. To avoid
  such points in the analysis of the light curves, we calculated the median of the photometric errors in each light curve and 
  removed those points that have photometric errors larger than ten times the median error.
  The final multi-band light curves of the objects are given in 
  Figures <ref>, <ref>, <ref> and <ref>.
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  § ANALYSIS
  
  
  
   §.§ Flux variability
  
  
  To characterise the flux variability nature of the sources, we calculated the intrinsic amplitude of 
  variability(σ_m) as outlined in <cit.>. The σ_m was 
  calculated from the observed light curves after subtracting the errors in the photometric measurements as
  follows (see also )
  
  
      Σ =√(1/n-1∑_i=1^N(m_i-<m>)^2),
  
  where <m> is the weighted mean of the magnitude in each band and σ_m is given by,
  
  
      σ_m =
        √(Σ^2 - ϵ^2)   if Σ≥ϵ,
      
               0,    otherwise
  
  Σ^2 is the variance of the light curve and ϵ represents the errors in the
  photometric measurements and it is calculated directly from individual measurement  errors ϵ_i as,
  
  
      ϵ^2 = 1/N∑_i=1^Nϵ_i^2
  
  
  The values of σ_m calculated for the sources along with the number of 
  data points available in each band are given in Table <ref>.  Considering 
  individual objects, σ_m in the infrared J and K bands is larger
  than the optical B-band for about 80% of the objects.  The mean  
  σ_m values along with their standard deviations in different bands for 
  FSRQ and BL Lac as a population are given in Table <ref>. 
  From Table <ref>, there are indications that the average σ_m in 
  the infrared J and K bands for the blazars studied here are larger than the 
  optical bands for both the FSRQ and BL Lac population, with the trend of increasing σ_m 
  towards longer wavelengths. However, given the large standard deviations to the 
  mean σ_m given in Table <ref>, it is difficult to make 
  a strong claim of σ_m increasing from shorter towards longer wavelengths
  for both FSRQ and BL Lac population. To ascertain the difference in the mean 
  σ_m between different wavelengths, we did a two sample Kolmogorov-Smironv (KS) test.  
  We restricted this test to B, R and J-bands, due to the availability of σ_m 
  measurements for all the sources in these bands.
  For BL Lacs, between B and R-bands, we found a D-statistics of 0.143 and a
  null hypothesis (no difference in σ_m between B and R bands) probability (p) 
  of 1.00.  Between B v/s J (R v/s J) we found D-value of 0.286(0.286) and null hypothesis
  p value of 0.883(0.883). Thus, for BL Lacs as a population, the observed hint of increased
  variability towards longer wavelengths is not statistically significant.
  Similarly, for FSRQs, between B and R-bands, our KS test gave a D-statistic of 0.167 and
  p of 0.760.  For variations between B and J-bands, we found a D-value of 0.400 and a 
  p of 0.011. And, between R and J-bands, the obtained D and p values are 0.233 and 0.342 respectively. Thus, the two sample KS test conducted between the distributions
  of σ_m between B and J-bands in FSRQs confirms that the two distributions are indeed different, and therefore, the difference in mean σ_m between J and B-band is
  statistically significant. Also, the differences in the mean σ_m between (B,R) and (R,J) are statistically insignificant.  Similarly in all the bands (except the K-band) the mean σ_m in BL Lacs is larger than that of FSRQs (see Table <ref>). Here too, as the standard deviations are larger, we did a KS test. For B-band, we found a D-statistics and p-value (null hypothesis probability of no difference in the σ_m between BL Lacs and FSRQs) of 0.324 and 0.504 respectively. Similarly, for V-band (R-band) we found D and p values of 0.412 (0.224) and 0.519 (0.901) respectively. For J-band we found a D and p of 0.314 and 0.543, while for K-band the D and p values are 0.289 and 0.834 respectively. Thus, from KS test we found that the distribution of σ_m between FSRQs and BL Lacs in all the bands is indistinguishable.  Analysis similar to the one carried out in this work is needed on a larger number of FSRQs and BL Lacs to confirm if BL Lacs indeed show σ_m larger than FSRQs on month to year like time scales.
  
  
   
  
  
  
  
  
  
  
  
  
   §.§ Colour Variability
  
  In addition to flux variations, blazars are also known to show spectral variations. The photometric measurements
  in BVRJK bands analysed here have contributions from both thermal emission coming from the accretion disk  and non-thermal emission coming from the relativistic jet. Both these processes contribute differently to the emission in different bands. Thus by studying spectral variations, one can distinguish the different components contributing to the observed flux variations. To characterise the nature of spectral variations in our sample of sources, we generated colour (B-J) -  magnitude (J-band) diagrams. Linear
  least squares fit was carried out on the data points in the colour-magnitude diagram by properly taking into account the errors in both colour and magnitude. For this we used BCES method [https://github.com/rsnemmen/BCES] <cit.> of orthogonal least squares to find the slope and intercept. A positive slope in the colour-magnitude diagram is an indication of a BWB 
  trend, while a negative slope is an indication of a RWB trend. The generated colour-magnitude diagrams for all the objects in our sample, along with the linear least squares fit are given in Figures
  <ref>, <ref> and <ref>. The results of the linear least squares fit to the 
  colour magnitude diagram are given in Table <ref>. 
  We considered BWB color-magnitude variation in 
  a source if the Spearman rank correlation coefficient is positive and the probability p of no correlation is lesser than 0.05. Similarly, a source is considered to have shown a RWB trend if values of the Spearman rank correlation coefficient is negative and p is lesser than 0.05. Using the above criteria we found that a large fraction of the sources in our sample ( 80%) showed a RWB trend. No spectral change with brightness was noticed in three sources, while another two sources showed a BWB trend. In two sources namely A0 0235+164 and 3C 379, we noticed both BWB and RWB trend. They showed a RWB trend upto a certain J-band brightness, and for magnitudes brighter than that they showed a BWB behaviour. For two sources PKS 0637-75 and BZQ J0850-1213 two deviant point each (shown as purple points in <ref> and <ref> respectively) were excluded from the fit, as inclusion of those points gave large errors in the derived slope and intercept. These observations indicate that complex spectral variations
  with brightness are invariably present in FSRQs and BL Lacs.
  
   
            
  
  
  
  
  
  
   §.§ Multi-wavelength cross-correlation
     
  All the sources analysed in this work are highly variable in both optical and infrared bands. 
  To check for correlation between the flux variations in different bands we used 
  the discrete correlation function (DCF; ) and the interpolated cross-correlation function (ICCF) 
  techniques <cit.>.
  The DCF method calculates the correlation function from unevenly sampled data. For two light 
  curves, say, x_i and y_j this method calculates a set of unbinned discrete correlation functions (UDCFs) given as
  
      UDCF_ij = (x_i - x) (y_j - y)/√((σ_x^2 - e_x^2)(σ_y^2 - e_y^2))
  
  here x and y are the mean of the two light curves, σ_x and σ_y are
  the standard deviations of the two light curves and e_x and e_y are the errors in each data point. For 
  any time lag (τ), DCF  is taken as the average of all UDCFs having the same τ. Taking the average of 
  M pairs of UDCFs for which τ - Δτ/2 ≤Δτ_ij < τ + Δτ/2, we get
  
      DCF(τ) = 1/M∑ UDCF_ij
  
  The error in each DCF bin is given as
  
      σ(τ) = 1/M-1∑[UDCF_ij - DCF(τ)]^2^1/2
  
  
  To estimate the lag we used 
  the centroid (τ_cen) of the cross-correlation function (CCF) defined as
  
  
      τ_cen = ∑_i τ_i CCF_i/∑_i CCF_i
  
  
   For calculating the centroid, we included those points in the CCF that are within 80% of the CCF peak. To get the errors in the CCF, we followed the Monte Carlo approach based on flux randomization (FR) and random subset selection (RSS) outlined in <cit.> that takes into account the time sampling of the light curves and the randomness in the flux measurements. For each pair of light curves using both FR and RSS, we carried out 10,000 iterations, and for each iteration we determined the τ_cen. Using these values, we generated the distribution of CCF values. The median of the CCF distribution is taken as the value of the lag. As the distribution of τ_cen obtained from 10,000 iterations is non-Gaussian, the upper (Δτ_up) and lower (Δτ_low)  uncertainties in τ_cen were determined such that 15.87% of iterations have τ > τ_cen + Δτ_up and 15.87% of iterations have τ < τ_cen - Δτ_low. The upper and lower error thus obtained in τ_cen  corresponds to ± 1 σ for a function that has a Gaussian distribution. The above procedures to obtain the lag and the associated error using DCF method was also repeated  using ICCF technique, with the exception that the light  curves were first brought to uniform time binning using linear interpolation. Thus, for all the sources we have two measurements of lag, one using the DCF method and the other one using ICCF method. The cross-correlation plots obtained between B and J bands is shown in Figures <ref> and <ref>. For majority of the sources, we found no lag between the flux variations in B and J-bands. Only for three sources, namely, PKS 1144-379 (lag = 19.3^+8.9_-10.7 days), PKS B1424-418 (lag= 13.4^+2.0_-3.4 days) and 3C 273 (lag = 461.0^+37.6_-30.2 days), we found statistically significant lags between flux variations in B and J-bands. 
  
  
  
  
  
  
  
  
  § DISCUSSION
  
  
  
   §.§ Flux variability
  
  The long term optical/UV and infrared flux variations observed in blazars could be due to intrinsic processes, extrinsic processes or a combination of both. Internal
  physical process include  those related to accretion disk caused by change in accretion rate <cit.> and those related to the relativistic jet such as blobs of emission moving through the helical magnetic fields of the jet <cit.>, the propagation of 
  shock waves in helical jets <cit.> or changes in the direction and speed of the jet flow <cit.>. Alternatively, the observed long term flux changes in blazars could also be caused by extrinsic geometrical effects associated with a change in
  our viewing angle to the moving blob of emission, that leads to variable Doppler boosting of the emitted radiation <cit.>. 
  The observed brightness in the UV, optical and the infrared bands in blazars is a combination of blue thermal emission from the accretion disk and red relativistically beamed non-thermal emission from their relativistic jets
  <cit.>. The 
  emission from the disk peaks in the UV band (the big blue bump at around 
  3000 Å) and its contribution decreases towards optical and infrared 
  wavelengths. In the infrared band, the thermal emission is thought to come
  from the dusty molecular torus that reprocesses the photons from the broad
  line region and the accretion disk into the infrared region <cit.>. As the thermal emission comes from a larger region in the accretion disk, they are expected to vary less 
  than the variations caused from the jet. The blazars studied here are γ-ray loud sources detected by the Fermi Gamma-ray Space Telescope
  <cit.>, and are believed to have strong relativistic jets whose non-thermal emission 
  dominates over the thermal accretion disk emission most of the time. Therefore, in this scenario, it is quite natural to expect σ_m to increase with wavelength. This
  is what is observed on individual sources in this work, where in majority of the objects (about 80%), σ_m in the infrared J and K bands is larger than that in the 
  optical B-band irrespective of FSRQs or BL Lacs. 
  Considering FSRQs and BL Lacs as a population,  the mean σ_m is found to progressively increase towards longer wavelengths in both FSRQs and BL Lacs 
  (see Table <ref>) though the standard deviations on each of the mean σ_m values are very large. A two sample KS test was therefore carried out on the distribution of σ_m values between different sets of filters such as (B,R), (R,J) and (B,J)
  for both FSRQs and BL Lacs. Statistically significant difference was only found between the 
  distribution of σ_m values in B and J-bands for FSRQs. 
  
  
  
  
  Among blazars, a prominent accretion disk 
  component is mostly seen in FSRQs that have broad emission lines and a radiatively efficient accretion disk relative to BL Lacs that have weak emission lines and radiatively inefficient accretion disk 
  <cit.>. However, it is not always true that accretion disk emission is prominent in FSRQs, 
  as there are instances wherein at increased activity states, the jet emission dominates the accretion disk emission making the accretion disk emission 
  non-distinguishable in the broad band SED. 
  Similarly, in BL Lacs too there are instances where the strong disc component is prominent in 
  their broad band SED <cit.> particularly during their faint activity states. 
  Irrespective of the conspicuous presence of accretion disk component in the broad band
  SEDs of both FSRQs and BL Lacs, the observed optical and infrared emission in most of our sample
  of γ-ray loud blazars studied here is unambiguously dominated by the red jet emission.
  Therefore, the increased σ_m towards longer wavelengths seen in majority of the 
  blazars studied here is understood by the dominance of the jet emission in them.
  Dividing the 37 blazars studied here based on the position of the synchrotron peak in their 
  broad band SED, we have 33 low synchrotron peaked sources (LSPs), 2 intermediate synchrotron peaked sources (ISPs) and 2 high synchrotron peaked sources (HSPs). In both the ISPs and HSPs too, the amplitude of variability in the infrared J-band is found to be larger than the optical B-band. Considering the average variability properties of the sample of FSRQs and BL Lacs as a whole, we found hints of increased variability amplitude towards longer 
  wavelengths in both the populations which again supports the dominance of jet emission in both the population of blazars.  
  The number of BL Lacs used in this study is only seven, therefore, similar analysis on a large number of FSRQs and
  BL Lacs are needed  to statistically confirm the indications observed from the mean σ_m values shown in Table <ref>. Considering, individual sources, 
  in a small fraction of the objects in the sample (that includes FSRQs and BL Lacs) the amplitude of variability is found to increase towards shorter wavelengths (see Table <ref>) which could be due to the low jet activity of the sources and the emission from the accretion
  disk becoming prominent over the jet emission. 
  
  Based on an analysis of a small sample of blazars utilising data from SMARTS, but spanning only two years <cit.> found FSRQs to show increased variability towards longer wavelengths, while they found BL Lacs not to show larger variability amplitude towards longer wavelengths. Similarly, <cit.> on analysis of a large number of FSRQs and BL Lacs found FSRQs to show larger σ_m in 
  the infrared J-band relative to the optical R-band. In the case of BL Lacs they found the σ_m to be larger in the R-band relative to J-band. However, in this work, from an analysis of about 10 years of data, in a larger number of FSRQs and BL Lacs we found hints of  increased σ_m towards longer wavelengths in both FSRQs and BL Lac, though statistically significant difference
  was only found between B and J-bands in FSRQs.  
  Our results for FSRQs are in general in agreement with that of <cit.> and <cit.>, however our
   results on BL Lacs are at odds with what is available in 
  literature <cit.>.
  
  Considering the sample of FSRQs and BL Lacs studied in this work as a whole, we found BL Lacs to have  larger mean σ_m in all the bands analysed here except K-band, though the standard deviations in the mean σ_m values are larger. Given the large standard deviation, it is difficult to conclusively say
  on strong statistical ground that BL Lacs show larger σ_m than FSRQs. Given this caveat, it is likely that both FSRQs and BL Lacs have the same intrinsic jet variability. However, as BL Lacs have small inclination than FSRQs, the difference in the σ_m values of optical and infrared variability  between FSRQs and BL Lacs with BL Lacs showing more σ_m compared
  to FSRQs could be due to BL Lacs viewed at smaller angles than FSRQs and subsequently subjected to higher
  Doppler boosting. Therefore, the difference in the mean σ_m between FSRQs and BL Lacs found here
  is in line with what is expected in the Unification scheme of radio-loud AGN, but to place
  the results that are found to be indicative here on a firm footing, this analysis needs
  to be extended on a large sample of FSRQs and BL Lacs.
  
  
  
   §.§ Colour Variability
  
  The observed optical, infrared emission in blazars is the sum of the blue thermal emission from the accretion disk and red non-thermal emission from the relativistic jet. The variations seen in the optical and infrared emission from blazars are often accompanied by colour(spectral) variations. Studies of the correlation between colour and brightness variations in blazars can give clues to the relation between accretion disk and jet emission in blazars.
  However the correlation between colour and brightness in blazars is highly debated in literature. Reports available  in literature indicate that among blazars, in BL Lacs the colour becomes bluer with increasing brightness (BWB;), 
  while in FSRQs the colour becomes redder with increasing brightness (RWB; ). Such a distinct colour magnitude relation between FSRQs and BL Lacs lead to speculate that the jets of FSRQs and BL Lacs are fundamentally different.
  
  However, in addition to the distinct colour magnitude relation between FSRQs and BL Lacs, varied colour magnitude relations too were seen in blazars. For example in the FSRQ
  3C 454.3, BWB trend was observed at certain epochs, while RWB trend was observed at certain other epochs and at one epoch both BWB and RWB trends were observed <cit.>. Also, in another FSRQ, 3C 345, BWB and RWB trends were simultaneously observed 
  <cit.>. In the BL Lac object  S5 0716+714, <cit.> 
  found  the colour to be weakly correlated with brightness with one clear BWB trend at a certain period, while <cit.> found BWB trend on both 
  inter-night and intra-night timescales. A RWB trend arises if the relativistic jet emission dominates over the less 
  variable thermal component from the accretion disk <cit.>.
  Alternatively, a BWB trend is attributed to increased amplitude of variability at shorter wavelengths
  <cit.>. Such a scenario, in the one zone synchrotron emission model could happen due to the injection of fresh electrons that have an energy distribution harder than that of the earlier cooler elections, that causes an increase in flux with
  a BWB trend <cit.>. Another explanation for a BWB behaviour is attributed to a change in the Doppler factor <cit.>. 
  According to <cit.> the BWB colour variation could be due to changes in the 
  Doppler factor, that could happen due to geometrical effects on the variation in the viewing angle of a curved and inhomogeneous jet.
  
  Majority of the sources in our sample showed a RWB trend on analysis of their complete data set. In all such sources, the amplitude of variability in the infrared J-band
  is larger than the optical B-band except three sources namely PKS 0537-441, PKS 0637-75 and 3C 273.
  Of these three, PKS 0537-441 is a BL Lac object while the other two are FSRQs. In two sources in our sample, namely A0 0235+164 and 3C 279 we found  a RWB trend upto a certain J-band brightness beyond which the spectrum
  changes shape with increasing brightness. Based on the position of the synchrotron peak, of the 37 blazars, 33 are LSP sources, 2 are HSPs and 2 are ISPs. The RWB trend is invariably found among all 
  the blazar sub-classes as well as among BL Lacs and FSRQs.
  All the objects studied here are γ-ray emitters and thus have conspicuous 
  jets in them. In this work we found similar flux and colour variability patterns in both FSRQs and BL Lacs. 
  This lead us to conclude that the jets in both FSRQs and BL Lacs have a similar 
  flux and colour variability behaviour.
  
  
  § SUMMARY
  
  
  We carried out an analysis of the multiwavelength (BVRJK) flux variability 
  characteristics of a sample of 37 blazars (that 
  includes 30 FSRQs and 7 BL Lacs)  taken
  from the SMARTS archives using data that spans over 10 years from 2008-2018. The main findings are summarized below.
  
  
  
    * All the sources analysed here showed long term flux variations in both the optical and infrared bands. In a majority of the sources (∼80%), the amplitude of variations in the infrared J-band is larger than the amplitude of variations in the optical B-band. Considering FSRQs and BL Lacs as a population, there 
  are hints that the  amplitude of flux variations in the longer
  wavelengths is larger than that at the shorter wavelengths. This behaviour is also
  evident among the different sub-classes of blazars, namely LSP, ISP and HSP sources. 
  For FSRQs we found statistically significant difference between the flux variations in 
  B and J bands. Such increased variations at the longer wavelengths is due to the dominance of 
  the jet emission in blazars.
  
  
    * Between FSRQs and BL Lacs as a population, the data analysed in this work is 
  suggestive (though statistically insignificant) of BL Lacs showing increased amplitude of flux variations in the optical 
  and infrared bands than FSRQs. This is understood as the jets of BL Lacs are 
  viewed at much smaller angles than FSRQs and thereby
  subjected to large variability amplitude due to Doppler boosting.
  
  
    * We found close correlations between colour and brightness changes with 
  a majority of FSRQs and BL Lacs predominantly showing a RWB trend.
  
  
    * The observed flux and colour variations in this work, lead to conclude 
  that the relativistic jets in FSRQs and BL Lacs are similar.
  
   
  
  
  
  § ACKNOWLEDGEMENTS
  
  We thank the referee for his/her critical comments that helped to improve the manuscript. This paper has made use of up-to-date SMARTS optical/near-infrared light curves that are available at www.astro.yale.edu/smarts/glast/home.php. Safna, P.Z., thanks CHRIST (Deemed to be University), for providing the facilities to carry out this work. The help by Arun Roy, Ujjwal Krishnan, Sabik PS, Raghu Krishnan and Shakhil PG is also thankfully acknowledged. 
  
  
  
  § DATA AVAILABILITY
  
  The photometric data underlying this article are publicly available from the Small and Moderate
  Aperture Research Telescope System (SMARTS[http://www.astro.yale.edu/smarts/fermi]) portal. 
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  @urlcharsothermakeother $&#_%
  @doi@urlcharsother ifnextchar [ @doi@
    @doi@[]
  @doi@[#1]#2tempa#1tempaempty http://dx.doi.org/#2 doi:#2http://dx.doi.org/#2 #1
  @eprint#1#2@eprint@#1:#2::nil
  @eprint@arXiv#1http://arxiv.org/abs/#1 arXiv:#1
  @eprint@dblp#1http://dblp.uni-trier.de/rec/bibtex/#1.xml
    dblp:#1
  @eprint@#1:#2:#3:#4niltempa #1tempb #2tempc
    #3tempc empty tempc tempb tempb tempa tempb empty tempb arXivifundefined
    mn@eprint@tempbtempb:tempcmn@eprint@tempbtempc
  
  [Ackermann et al.,Ackermann
    et al.2015]2015ApJ...810...14A
  Ackermann M.,  et al., 2015, @doi [] 10.1088/0004-637X/810/1/14,
    https://ui.adsabs.harvard.edu/abs/2015ApJ...810...14A 810, 14
  
  [Angel & StockmanAngel &
    Stockman1980]1980ARA   A..18..321A
  Angel J. R. P.,  Stockman H. S.,  1980, @doi []
    10.1146/annurev.aa.18.090180.001541, https://ui.adsabs.harvard.edu/abs/1980ARA   A..18..321A 18, 321
  
  [AntonucciAntonucci1984]Antonucci
  Antonucci R. R. J.,  1984, @doi [] 10.1086/161816, http://adsabs.harvard.edu/abs/1984ApJ...278..499A 278, 499
  
  [AntonucciAntonucci1993]1993ARA   A..31..473A
  Antonucci R.,  1993, @doi [] 10.1146/annurev.aa.31.090193.002353,
    https://ui.adsabs.harvard.edu/abs/1993ARA   A..31..473A 31, 473
  
  [Bonning et al.,Bonning
    et al.2012]Bonning2012ApJ...756...13B
  Bonning E.,  et al., 2012, @doi [] 10.1088/0004-637X/756/1/13, http://adsabs.harvard.edu/abs/2012ApJ...756...13B 756, 13
  
  [BöttcherBöttcher2007]2007Ap   SS.309...95B
  Böttcher M.,  2007, @doi [] 10.1007/s10509-007-9404-0, https://ui.adsabs.harvard.edu/abs/2007Ap   SS.309...95B 309, 95
  
  [Camenzind & KrockenbergerCamenzind &
    Krockenberger1992]1992A   A...255...59C
  Camenzind M.,  Krockenberger M.,  1992, , https://ui.adsabs.harvard.edu/abs/1992A   A...255...59C 255, 59
  
  [D'Amicis, Nesci, Massaro, Maesano,
    Montagni  & D'AlessioD'Amicis et al.2002]2002PASA...19..111D
  D'Amicis R.,  Nesci R.,  Massaro E.,  Maesano M.,  Montagni F.,
    D'Alessio F.,  2002, @doi [] 10.1071/AS01113, https://ui.adsabs.harvard.edu/abs/2002PASA...19..111D 19, 111
  
  [Dermer & SchlickeiserDermer &
    Schlickeiser1993]1993ApJ...416..458D
  Dermer C. D.,  Schlickeiser R.,  1993, @doi [] 10.1086/173251,
    https://ui.adsabs.harvard.edu/abs/1993ApJ...416..458D 416, 458
  
  [Edelson & KrolikEdelson &
    Krolik1988]Edilson1988ApJ...333..646E
  Edelson R. A.,  Krolik J. H.,  1988, @doi [] 10.1086/166773, http://adsabs.harvard.edu/abs/1988ApJ...333..646E 333, 646
  
  [Fiorucci, Ciprini  & TostiFiorucci
    et al.2004]2004A   A...419...25F
  Fiorucci M.,  Ciprini S.,   Tosti G.,  2004, @doi []
    10.1051/0004-6361:20034218, https://ui.adsabs.harvard.edu/abs/2004A   A...419...25F 419, 25
  
  [Gaskell & PetersonGaskell &
    Peterson1987]1987ApJS...65....1G
  Gaskell C. M.,  Peterson B. M.,  1987, @doi [] 10.1086/191216,
    https://ui.adsabs.harvard.edu/abs/1987ApJS...65....1G 65, 1
  
  [Gaskell & SparkeGaskell &
    Sparke1986]1986ApJ...305..175G
  Gaskell C. M.,  Sparke L. S.,  1986, @doi [] 10.1086/164238, https://ui.adsabs.harvard.edu/abs/1986ApJ...305..175G 305, 175
  
  [Gaur et al.,Gaur
    et al.2019]2019MNRAS.484.5633G
  Gaur H.,  et al., 2019, @doi [] 10.1093/mnras/stz322, https://ui.adsabs.harvard.edu/abs/2019MNRAS.484.5633G 484, 5633
  
  [GhiselliniGhisellini2019]2019arXiv191111777G
  Ghisellini G.,  2019, arXiv e-prints, https://ui.adsabs.harvard.edu/abs/2019arXiv191111777G p. arXiv:1911.11777
  
  [Ghisellini & CelottiGhisellini &
    Celotti2001]2001A   A...379L...1G
  Ghisellini G.,  Celotti A.,  2001, @doi []
    10.1051/0004-6361:20011338, https://ui.adsabs.harvard.edu/abs/2001A   A...379L...1G 379, L1
  
  [Ghisellini, Maraschi  &
    TavecchioGhisellini et al.2009]2009MNRAS.396L.105G
  Ghisellini G.,  Maraschi L.,   Tavecchio F.,  2009, @doi []
    10.1111/j.1745-3933.2009.00673.x, https://ui.adsabs.harvard.edu/abs/2009MNRAS.396L.105G 396, L105
  
  [Gopal-Krishna & WiitaGopal-Krishna &
    Wiita1992]1992A   A...259..109G
  Gopal-Krishna Wiita P. J.,  1992, , https://ui.adsabs.harvard.edu/abs/1992A   A...259..109G 259, 109
  
  [Jones, O'dell  & SteinJones
    et al.1974]1974ApJ...188..353J
  Jones T. W.,  O'dell S. L.,   Stein W. A.,  1974, @doi []
    10.1086/152724, https://ui.adsabs.harvard.edu/abs/1974ApJ...188..353J 188, 353
  
  [Kirk, Rieger  & MastichiadisKirk
    et al.1998]1998A   A...333..452K
  Kirk J. G.,  Rieger F. M.,   Mastichiadis A.,  1998, , https://ui.adsabs.harvard.edu/abs/1998A   A...333..452K 333, 452
  
  [KoniglKonigl1981]Konigl1981ApJ...243..700K
  Konigl A.,  1981, @doi [] 10.1086/158638, http://adsabs.harvard.edu/abs/1981ApJ...243..700K 243, 700
  
  [Larionov, Villata  &
    RaiteriLarionov et al.2010]2010A   A...510A..93L
  Larionov V. M.,  Villata M.,   Raiteri C. M.,  2010, @doi []
    10.1051/0004-6361/200913536, https://ui.adsabs.harvard.edu/abs/2010A   A...510A..93L 510, A93
  
  [Larionov et al.,Larionov
    et al.2013]2013ApJ...768...40L
  Larionov V. M.,  et al., 2013, @doi [] 10.1088/0004-637X/768/1/40,
    https://ui.adsabs.harvard.edu/abs/2013ApJ...768...40L 768, 40
  
  [Li & CaoLi &
    Cao2008]2008MNRAS.387L..41L
  Li S.-L.,  Cao X.,  2008, @doi []
    10.1111/j.1745-3933.2008.00480.x, https://ui.adsabs.harvard.edu/abs/2008MNRAS.387L..41L 387, L41
  
  [Liodakis, Romani, Filippenko,
    Kiehlmann, Max-Moerbeck, Readhead  & ZhengLiodakis
    et al.2018]2018MNRAS.480.5517L
  Liodakis I.,  Romani R. W.,  Filippenko A. V.,  Kiehlmann S.,
    Max-Moerbeck W.,  Readhead A. C. S.,   Zheng W.,  2018, @doi
    [] 10.1093/mnras/sty2264, https://ui.adsabs.harvard.edu/abs/2018MNRAS.480.5517L 480, 5517
  
  [Lynden-BellLynden-Bell1969]Lynden1969Natur.223..690L
  Lynden-Bell D.,  1969, @doi [] 10.1038/223690a0, http://adsabs.harvard.edu/abs/1969Natur.223..690L 223, 690
  
  [Marscher et al.,Marscher
    et al.2008]2008Natur.452..966M
  Marscher A. P.,  et al., 2008, @doi [] 10.1038/nature06895, https://ui.adsabs.harvard.edu/abs/2008Natur.452..966M 452, 966
  
  [Mastichiadis & KirkMastichiadis &
    Kirk2002]2002PASA...19..138M
  Mastichiadis A.,  Kirk J. G.,  2002, @doi [] 10.1071/AS01108,
    https://ui.adsabs.harvard.edu/abs/2002PASA...19..138M 19, 138
  
  [Meng, Zhang, Wu, Ma  & ZhouMeng
    et al.2018]2018ApJS..237...30M
  Meng N.,  Zhang X.,  Wu J.,  Ma J.,   Zhou X.,  2018, @doi []
    10.3847/1538-4365/aacffe, https://ui.adsabs.harvard.edu/abs/2018ApJS..237...30M 237, 30
  
  [Mücke & ProtheroeMücke &
    Protheroe2001]2001APh....15..121M
  Mücke A.,  Protheroe R. J.,  2001, @doi [Astroparticle Physics]
    10.1016/S0927-6505(00)00141-9, https://ui.adsabs.harvard.edu/abs/2001APh....15..121M 15, 121
  
  [Mücke, Protheroe, Engel, Rachen
    & StanevMücke et al.2003]2003APh....18..593M
  Mücke A.,  Protheroe R. J.,  Engel R.,  Rachen J. P.,   Stanev
    T.,  2003, @doi [Astroparticle Physics] 10.1016/S0927-6505(02)00185-8,
    https://ui.adsabs.harvard.edu/abs/2003APh....18..593M 18, 593
  
  [Papadakis, Villata  &
    RaiteriPapadakis et al.2007]2007A   A...470..857P
  Papadakis I. E.,  Villata M.,   Raiteri C. M.,  2007, @doi []
    10.1051/0004-6361:20077516, https://ui.adsabs.harvard.edu/abs/2007A   A...470..857P 470, 857
  
  [Perlman, Addison, Georganopoulos,
    Wingert  & GraffPerlman et al.2008]2008bves.confE...9P
  Perlman E.,  Addison B.,  Georganopoulos M.,  Wingert B.,   Graff P.,
     2008, in Blazar Variability across the Electromagnetic Spectrum. p. 9
    (@eprint arXiv 0807.2119)
  
  [Peterson, Wanders, Horne, Collier,
    Alexander, Kaspi  & MaozPeterson
    et al.1998]petersonccfA1998PASP..110..660P
  Peterson B. M.,  Wanders I.,  Horne K.,  Collier S.,  Alexander T.,
    Kaspi S.,   Maoz D.,  1998, @doi [] 10.1086/316177, http://adsabs.harvard.edu/abs/1998PASP..110..660P 110, 660
  
  [Pollack, Pauls  & WiitaPollack
    et al.2016]2016ApJ...820...12P
  Pollack M.,  Pauls D.,   Wiita P. J.,  2016, @doi []
    10.3847/0004-637X/820/1/12, https://ui.adsabs.harvard.edu/abs/2016ApJ...820...12P 820, 12
  
  [Raiteri et al.,Raiteri
    et al.2003]2003A   A...402..151R
  Raiteri C. M.,  et al., 2003, @doi [] 10.1051/0004-6361:20030256,
    https://ui.adsabs.harvard.edu/abs/2003A   A...402..151R 402, 151
  
  [Raiteri et al.,Raiteri
    et al.2008]2008A   A...480..339R
  Raiteri C. M.,  et al., 2008, @doi [] 10.1051/0004-6361:20079044,
    https://ui.adsabs.harvard.edu/abs/2008A   A...480..339R 480, 339
  
  [Rajput, Stalin, Sahayanathan, Rakshit
    & MandalRajput et al.2019]2019MNRAS.486.1781R
  Rajput B.,  Stalin C. S.,  Sahayanathan S.,  Rakshit S.,   Mandal
    A. K.,  2019, @doi [] 10.1093/mnras/stz941, https://ui.adsabs.harvard.edu/abs/2019MNRAS.486.1781R 486, 1781
  
  [Rakshit & StalinRakshit &
    Stalin2017]Rakshit2017ApJ...842...96R
  Rakshit S.,  Stalin C. S.,  2017, @doi []
    10.3847/1538-4357/aa72f4, http://adsabs.harvard.edu/abs/2017ApJ...842...96R 842, 96
  
  [Rakshit, Stalin, Muneer, Neha  &
    PaliyaRakshit et al.2017]2017ApJ...835..275R
  Rakshit S.,  Stalin C. S.,  Muneer S.,  Neha S.,   Paliya V. S.,
    2017, @doi [] 10.3847/1538-4357/835/2/275, https://ui.adsabs.harvard.edu/abs/2017ApJ...835..275R 835, 275
  
  [Rani, Stalin  & RakshitRani
    et al.2017]2017MNRAS.466.3309R
  Rani P.,  Stalin C. S.,   Rakshit S.,  2017, @doi []
    10.1093/mnras/stw3228, https://ui.adsabs.harvard.edu/abs/2017MNRAS.466.3309R 466, 3309
  
  [ReesRees1984]Rees1984ARA   A..22..471R
  Rees M. J.,  1984, @doi [] 10.1146/annurev.aa.22.090184.002351,
    http://adsabs.harvard.edu/abs/1984ARA26A..22..471R 22, 471
  
  [Sarkar et al.,Sarkar
    et al.2019]2019ApJ...887..185S
  Sarkar A.,  et al., 2019, @doi [] 10.3847/1538-4357/ab5281, https://ui.adsabs.harvard.edu/abs/2019ApJ...887..185S 887, 185
  
  [Sesar et al.,Sesar
    et al.2007]2007AJ....134.2236S
  Sesar B.,  et al., 2007, @doi [] 10.1086/521819, https://ui.adsabs.harvard.edu/abs/2007AJ....134.2236S 134, 2236
  
  [Sikora, Begelman  & ReesSikora
    et al.1994]1994ApJ...421..153S
  Sikora M.,  Begelman M. C.,   Rees M. J.,  1994, @doi []
    10.1086/173633, https://ui.adsabs.harvard.edu/abs/1994ApJ...421..153S 421, 153
  
  [Stalin et al.,Stalin
    et al.2009]2009MNRAS.399.1357S
  Stalin C. S.,  et al., 2009, @doi []
    10.1111/j.1365-2966.2009.15354.x, https://ui.adsabs.harvard.edu/abs/2009MNRAS.399.1357S 399, 1357
  
  [Ulrich, Maraschi  & UrryUlrich
    et al.1997]Ulrich1997ARA   A..35..445U
  Ulrich M.-H.,  Maraschi L.,   Urry C. M.,  1997, @doi []
    10.1146/annurev.astro.35.1.445, http://adsabs.harvard.edu/abs/1997ARA26A..35..445U 35, 445
  
  [Urry & MushotzkyUrry &
    Mushotzky1982]UrryMushosky1982ApJ...253...38U
  Urry C. M.,  Mushotzky R. F.,  1982, @doi [] 10.1086/159607, http://adsabs.harvard.edu/abs/1982ApJ...253...38U 253, 38
  
  [Urry & PadovaniUrry &
    Padovani1995]Urry1995PASP..107..803U
  Urry C. M.,  Padovani P.,  1995, @doi [] 10.1086/133630, http://adsabs.harvard.edu/abs/1995PASP..107..803U 107, 803
  
  [Vagnetti, Trevese  & NesciVagnetti
    et al.2003]2003ApJ...590..123V
  Vagnetti F.,  Trevese D.,   Nesci R.,  2003, @doi []
    10.1086/374889, https://ui.adsabs.harvard.edu/abs/2003ApJ...590..123V 590, 123
  
  [Villata et al.,Villata
    et al.2004]2004A   A...421..103V
  Villata M.,  et al., 2004, @doi [] 10.1051/0004-6361:20035895, https://ui.adsabs.harvard.edu/abs/2004A   A...421..103V 421, 103
  
  [Villata et al.,Villata
    et al.2009]2009A   A...501..455V
  Villata M.,  et al., 2009, @doi [] 10.1051/0004-6361/200912065,
    https://ui.adsabs.harvard.edu/abs/2009A   A...501..455V 501, 455
  
  [Wagner & WitzelWagner &
    Witzel1995]Wagner1995ARA   A..33..163W
  Wagner S. J.,  Witzel A.,  1995, @doi []
    10.1146/annurev.aa.33.090195.001115, http://adsabs.harvard.edu/abs/1995ARA26A..33..163W 33, 163
  
  [White & PetersonWhite &
    Peterson1994]whiteprtrson1994PASP..106..879W
  White R. J.,  Peterson B. M.,  1994, @doi [] 10.1086/133456,
    http://adsabs.harvard.edu/abs/1994PASP..106..879W 106, 879
  
  [Wu, Zhou, Ma  & JiangWu
    et al.2011]2011MNRAS.418.1640W
  Wu J.,  Zhou X.,  Ma J.,   Jiang Z.,  2011, @doi []
    10.1111/j.1365-2966.2011.19565.x, https://ui.adsabs.harvard.edu/abs/2011MNRAS.418.1640W 418, 1640
  
  [Zhang, Zhou, Zhao  & DaiZhang
    et al.2015]2015RAA....15.1784Z
  Zhang B.-K.,  Zhou X.-S.,  Zhao X.-Y.,   Dai B.-Z.,  2015, @doi
    [Research in Astronomy and Astrophysics] 10.1088/1674-4527/15/11/002, https://ui.adsabs.harvard.edu/abs/2015RAA....15.1784Z 15, 1784
  
  [Akritas & Bershady1996]1996ApJ...470..706A Akritas M. G., Bershady M. A., 1996, ApJ, 470, 706
  
  [Nemmen et al.2012]2012Sci...338.1445N Nemmen R. S., Georganopoulos M., Guiriec S., Meyer E. T., Gehrels N., Sambruna R. M., 2012, Sci, 338, 1445
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  