Induced surface fluxes: A new framework for attributing

A new framework is presented for analysing the proximate causes of model Arctic sea ice biases, 7 demonstrated with the CMIP5 model HadGEM2-ES. In this framework the Arctic sea ice volume is treated as a 8 consequence of the integrated surface energy balance, via the volume balance. A simple model allows the local 9 dependence of the surface flux on specific model variables to be described, as a function of time and space. 10 When these are combined with reference datasets, it is possible to estimate the surface flux bias induced by the 11 model bias in each variable. The method allows the role of the surface albedo and ice thickness-growth 12 feedbacks in sea ice volume balance biases to be quantified, along with the roles of model bias in variables not 13 directly related to the sea ice volume. It shows biases in the HadGEM2-ES sea ice volume simulation to be due 14 to a bias in spring surface melt onset date, partly countered by a bias in winter downwelling longwave radiation. 15 The framework is applicable in principle to any model and has the potential to greatly improve understanding of 16 the reasons for ensemble spread in modelled sea ice state. A secondary finding is that observational uncertainty 17 is the largest cause of uncertainty in the induced surface flux bias calculation. 18


20
The Arctic sea ice cover has witnessed rapid change during the past 30 years, with a decline in September extent 21 of 1.05 x 10 6 km 2 / decade from 1986 to 2015 (HadISST1.2, Rayner et al 2003). In association with the changes in extent, Arctic sea ice thinning has been observed from submarine and satellite data (Rothrock et al, 2008, 23 Lindsay and Schweiger, 2015). Arctic sea ice has also become younger on average as older ice has been lost 24 (Maslanik et

32
Evaluating sea ice extent or volume with respect to reference datasets shows that some models reproduce 33 present-day sea ice state more accurately than others (e.g. Wang and Overland, 2012;Massonnet et al, 2012; does not necessarily imply an accurate future projection of sea ice change, as a correct simulation can be obtained by accident due to cancelling model errors. Sea ice extent in particular is known to be a very unsuitable While satellites have in many cases been able to produce Arctic-wide measurements of some characteristics, 23 most notably sea ice concentration, the relative lack of in situ observations against which these can be calibrated 24 means knowledge of the observational biases is limited. Reanalysis data over the Arctic is also more subject to 25 model errors than in other regions, due to errors in atmospheric forcing, and the existence of fewer direct 26 observations available for assimilation (Lindsay et al, 2014). The approach of this study is to use a wide range of 27 observational data to evaluate modelled sea ice state and surface radiative fluxes, setting results in the context of 28 in situ validation studies. The same datasets are then used as reference datasets for the induced surface flux 29 framework.

30
To evaluate modelled sea ice concentration, we use the HadISST1.2 dataset (Rayner et al, 2003), derived from 31 passive microwave observations. To evaluate modelled sea ice thickness Arctic-wide, we use the ice-ocean 32 model PIOMAS (Zhang and Rothrock, 2003), which is forced with the NCEP reanalysis and assimilates ice

31
This is likely to be associated with the ice thickness-growth feedback, whereby a steeper temperature gradient 32 induces stronger conduction and hence ice growth for thin ice. Indeed, negative correlations between summer 33 sea ice and sea ice growth the following winter are a ubiquitous feature of CMIP5 models (Massonnet et al,

35
The bias in ice volume balance is associated with a bias in ice energy uptake. In calculating this bias, and for the 36 rest of the study, we assume an ice density of 917 kgm -3 , the constant value used by HadGEM2-ES. For importance to the sea ice heat budget (Keen et al, 2018). Hence for the main purposes of this study, we 11 concentrate on the surface energy balance, and neglect the other two terms (although the ocean heat 12 convergence is briefly discussed in Section 6).

13
Surface radiative fluxes are now evaluated. In the following discussion, and throughout this study, the

25
As the June net SW bias is likely to result from a surface albedo bias we discuss the parameters affecting surface 26 albedo over sea ice in HadGEM2-ES: ice fraction, snow thickness and surface melt onset. Ice fraction has 27 already been evaluated and for snow thickness no reference dataset is available; however, surface melt onset can 28 be evaluated using satellite observations ( Figure 6). We define the date of melt onset for any grid cell as the first 29 day on which the surface temperature exceeds -1°C (varying this threshold by 0.5°C in either direction changes 30 the date in only a small minority of grid cells). The average date of melt onset as estimated by this method 31 ( Figure 5a) is then compared to that measured by the satellite-derived dataset described in Section 3 (Figure 5b), 32 with model bias shown in Figure 5c. Large spatial variability is evident in the observations. Melt onset occurs in 33 early to mid-May around the Arctic Ocean coasts, but much later in the Central Arctic, around mid-June. In 34 contrast, the HadGEM2-ES surface melt onset date is in mid-to late May across the Arctic Ocean, without the 35 strong gradients seen in the observations. This would cause a surface albedo, and hence net SW, bias with a strong maximum in the Central Arctic, similar to that discussed above.
downwelling LW, the mean model biases from December-April are -16, -22 and -40 Wm -2 for ERAI, CERES   1   and ISCCP-FD respectively; for upwelling LW, the biases are 11, 16 and 18 Wm -2 for CERES, ERAI    The surface radiation evaluation provides clues as to the causes of the HadGEM2-ES ice volume balance bias, 10 but also underlines why a more detailed analysis is required to properly quantify these causes. For example, in 11 winter the low downwelling LW bias provides a clear mechanism for the bias in ice freezing. However, the 12 reference datasets also suggest a low bias in upwelling LW that at first sight would tend to counteract this.. In 13 fact, these biases are fully consistent: a low bias in downwelling LW would be expected to cause a low surface 14 temperature bias, causing both a high bias in ice growth and a low bias in upwelling LW. The qualitative 15 evaluation fails to capture the full causal relationship; for this, an analysis of exactly how the downwelling LW bias affects surface flux, including the upwelling LW response, is necessary.
In the summer, meanwhile, the surface radiation evaluation suggests that a bias in net SW radiation is 18 responsible for the ice volume balance bias, and that this in turn is related to a surface albedo bias. However, at 19 least two possible drivers of this have been identified: the surface melt onset bias, and the underlying ice area 20 bias that is itself likely to be caused by the ice volume balance bias. Quantifying the extent to which each driver 21 is important, and at which times of year, is likely to help in resolving this circular causal loop. In the next section, it is described how the effects of each model bias on the sea ice volume balance can be separated and 23 quantified through their effect on the surface flux, in the ISF framework.  We are motivated by the observation that each of the model biases described above affects the sea ice volume 28 balance by acting through the total downwards surface energy flux (referred to as the surface flux). An excess of 29 downwelling radiation leads directly to a higher surface flux, higher sea ice energy uptake and a bias towards ice 30 melting. A bias in ice area, or in surface melt onset, is associated with a bias in surface albedo, hence a bias in 31 net SW, and in the total surface flux. Finally, a bias in ice thickness alters the thermodynamics of the entire 32 snow-ice column, and of the near-surface atmosphere layer. Flux continuity considerations show that in the Each of these relationships can be quantified, in principle, at any point in model space and time. Specifically, 1 the rate at which the surface flux depends on each variable alone, with others being held constant, can be 2 estimated. To this end, we approximate the surface flux sfc F at each point in model space

10
The function t x g , and its derivation are described fully in Appendix A, but are summarised briefly here. The   ice. The effect of estimating ISF bias at each point separately, and then averaging to determine large-scale 31 effects, is to bypass all nonlinearities.

32
A second advantage of this approach lies in the quasi-independence of the variables: while each variable may 33 affect the others over timescales varying from days to months, each affects the surface flux instantaneously (in t x g , , therefore, the sum of the ISF biases, over all variables, must approach the true model surface flux bias 1 (although it will be seen that this claim is impossible to evaluate precisely due to observational uncertainty). In 2 this way, large-scale model biases in surface flux, and hence sea ice volume balance, can be broken down into 3 separate contributions from model biases in each of the independent variables.

4
The ISF calculation process is now illustrated for two processes in turn. The model bias in melting surface 5 fraction for the month of June 1980 is positive over most of the Arctic, although only weakly so towards the 6 coasts ( Figure 6a), reflecting melt onset modelled earlier than observed during this month. The reference dataset 7 is derived from SSMI microwave observations. The rate of change of surface flux with respect to melt onset 8 occurrence tends to be higher in the Central Arctic (Figure 6b). This reflects a greater tendency to clear skies in 9 the Central Arctic, as this field is effectively downwelling SW multiplied by the difference in parameterised 10 albedos. The product of these two fields represents the modelled surface flux bias induced by the model bias in   Figure 7, together with total ISF bias, net radiative flux bias estimated by the direct 27 radiation evaluation relative to ISCCP-FD, CERES and ERAI, and also sea ice energy uptake biases implied by 28 the seasonal ice volume balance bias relative to PIOMAS. The ISF biases are also shown in Table 1   In both seasons, these results are consistent with the radiation and ice volume balance evaluation. equivalent to an extra 8cm of melt. This is associated with the overly fast retreat of sea ice in HadGEM2-ES, 1 and the low extents in late summer, as noted in Section 3. Thirdly, from October-March the downwelling LW season progresses, equivalent to total additional ice growth of 24cm.

6
Internal variability in the ISF biases is measured by taking the standard deviation of the whole-Arctic ISF bias 7 for each process and month across all 20 years in the model period, and all four ensemble members used.

8
Variability is highest in the ice area term, reaching 4.0 Wm -2 in July. Variability reaches considerable size in 9 some other terms in some months, for example 1.1 Wm -2 for surface melt onset in June, 1.9 Wm -2 for ice 10 thickness in November, but is otherwise mainly under 1 Wm -2 in magnitude. In each case, therefore, the ISF

26
The November-March ice thickness ISF bias displays a pattern which is almost identical to that of the August 27 ice fraction ISF bias, but with the opposite sign, with high negative values on the Atlantic side of the Arctic 28 rising to near-zero values in the Beaufort Sea. Hence this model bias has the reverse effect to that of the ice 29 concentration bias, reducing existing ice thickness biases by promoting additional ice growth in these areas.

30
Finally, the winter downwelling LW ISF bias is much more spatially uniform, but displays slightly higher 31 values on the Pacific side of the Arctic than the Atlantic side, a different pattern to that displayed by the 32 downwelling LW itself in Figure 4d. The contrast is due to the role the effective ice thickness scale factor plays 33 in determining the induced surface flux bias; the thicker ice present on the American side of the Arctic, tends to greatly reduce the flux bias. This represents the thickness-growth feedback; thicker ice will grow less quickly than thin ice under the same atmospheric conditions. The downwelling LW bias tends to increase ice growth 1 Arctic-wide, but less so in regions where ice is already thick.

2
The spatial patterns of total ISF bias shows many similarities to the total net radiation bias evaluated by CERES 3 in most months of the year (Figure 8), notably a tendency in July and August for positive surface flux biases to 4 be concentrated on the Atlantic side of the Arctic, and a tendency throughout the freezing season for negative 5 surface flux biases to be least pronounced in the Beaufort Sea, where the ice thickness biases are lowest. We 6 note that the spatial pattern of amplification of the ice thickness seasonal cycle displayed in Figure 3 is very 7 similar, with amplification most pronounced in the Atlantic sector, and least pronounced in the Beaufort Sea.

8
The surface flux biases produced by ice fraction biases in August, and ice thickness biases in November, 9 provide reasons for the spatial variation in amplification of the ice thickness seasonal cycle seen in Figure 4, as 10 well as the close resemblance of this pattern to the model biases in annual mean ice thickness. Ice which is 11 thinner in the annual mean will tend to melt faster in summer, due to the net SW biases associated with greater 12 creation of open water (the surface albedo feedback), and to freeze faster in winter, due to greater conduction of 13 energy through the ice (the ice thickness-growth feedback).

16
It is helpful to divide the processes examined into feedbacks (surface flux biases induced by biases in the sea ice 17 state itself) and forcings (those induced by variables external to the sea ice state). In this sense, a 'forcing' refers 18 to a variable which is independent of the sea ice volume on instantaneous timescales, rather than being used in

27
On the other hand, the ice thickness ISF bias (specifically during the freezing season) can be identified with the 28 effect of the thickness-growth feedback on the sea ice state. This is perhaps less obvious, as the ice thickness 29 affects the estimated surface flux via the surface temperature and upwelling LW radiation, while the thickness-30 growth feedback is usually understood to result from differences in conduction. However, the assumption of 31 flux continuity at the surface in constructing the estimated surface flux means that the cooler surface 32 temperatures, and shallower temperatures gradients occurring for thicker ice categories are manifestations of the 33 same process. Slower ice growth at higher ice thicknesses is caused by a smaller negative surface flux, and the 34 surface temperature is the mechanism by which this is demonstrated. Hence the effect of the thickness-growth   The ISF biases, summed over all independent variables, should approach the true total surface flux bias.

13
However, this is difficult to evaluate as the true surface flux bias is not known. Hence it is necessary to use 14 proxy quantities to evaluate the total ISF bias: directly evaluated surface net radiation bias (relative to ISCCP-15 FD, ERAI and CERES respectively); and ice energy uptake bias, derived from ice volume balance bias relative 16 to PIOMAS.

17
For most months of the year, all estimates of total ISF bias fall within the spread of these four datasets ( Figure  6), the exceptions being June and July when total ISF bias is smaller than all surface flux proxies. However, the 19 spread is extremely large. For example, in the month of January the estimates of total ISF bias are -12.3, -8.2 20 and -6.1 Wm -2 (with ISCCP-FD, CERES and ERAI used as downwelling radiation datasets respectively), while 21 the estimates of net radiation bias are -18.2, -11.6 and 0.6 Wm -2 from ISCCP-FD, CERES and ERAI 22 respectively, and ice heat uptake bias is estimated as -10.1 Wm -2 . Hence it is difficult to evaluate the total ISF 23 bias within current observational constraints, and at best it can be said that the total ISF bias is qualitatively 24 consistent, over the year as a whole, with the surface flux bias proxies. A possible cause of the lower total ISF 25 bias in June and July is the 'missing process' of snow on ice, which cannot be evaluated here due to the lack of a 26 reference dataset. The early surface melt onset, and sea ice fraction loss, as modelled by HadGEM2-ES, would 27 be associated with an early loss of snow on ice, with an additional surface albedo bias and hence an additional 28 ISF bias.

29
On the other hand, the annual mean ice heat uptake bias (0.0 Wm -2 ) provides a strong constraint on the annual 30 mean surface heat flux bias, in the absence of a significant oceanic heat convergence contribution. For example, 31 the annual mean total ISF biases are -3.6 Wm -2 and -4.5 Wm -2 when CERES and ISCCP-FD are used as 32 reference datasets respectively; these would imply sea ice thickening of 7m and 9m over the 1980-1999 in bias is related to the tendency for the ISF analysis to account for a greater bias towards ice growth in winter (45-multiple reference datasets whose errors are not constrained to correlate in a physically realistic sense, but may 1 also be related to the missing processes in June and July.

2
Observational error is one potential cause of error in the ISF biases. An idea as to the potential magnitude of this 3 can be seen from the large spread in SW and LW ISF bias (across different datasets) during summer ( Figure 6).

4
For example, in July the model downwelling LW bias with respect to ERAI produces an aggregated ISF bias of 5 -7.0 Wm -2 , but that with respect to CERES produces an aggregate ISF bias of 8.0 Wm -2 . Calculation of ice area 6 ISF biases using NSIDC and HadISST.2 as reference, described in section 2.2 and not shown here, showed a 7 similar magnitude of uncertainty in the ice area term (±10Wm -2 in summer, ±2Wm -2 in winter). We note that the 8 evidence from in situ validation studies suggests that the winter downwelling LW estimates of ERAI and 9 CERES are more likely to be accurate than that of ISCCP-FD. Hence the downwelling LW ISF bias is likely to 10 be estimated more accurately when ERAI or CERES are reference dataset, and the bias towards ice growth is 11 likely to lie closer to the lower end of the range (20cm).

12
In Appendix B, inherent theoretical errors in the ISF analysis are discussed and are found to be small relative to

33
(occurring mainly in areas of thinner ice, of greater importance in early winter) and by the downwelling LW bias (more spatially uniform, in late winter). It is unclear that any significant role is played by the downwelling the induced surface flux biasis more than balanced by that induced by downwelling LW. However this may 1 have a role in causing the later melt onset bias, as discussed below.

2
The means by which the external forcings -anomalous LW winter cooling, and early late spring melt onset -3 cause an amplified seasonal cycle in sea ice thickness are clear. It can also be seen how, in the absence of other 4 forcings, these combine to create an annual mean sea ice thickness that is biased low, as seen in Section 4. The 10 Acting together, the ice thickness-growth feedback and surface albedo feedback create a strong association 11 between lower ice thicknesses and amplified seasonal cycles, because ice which tends to be thinner will both 12 grow faster during the winter, and melt faster during the summer. Hence the melt onset bias, acting alone, would 13 induce a seasonal cycle of sea ice thickness lower in the annual mean, but also more amplified, than that

23
The feedbacks of the sea ice state explain the association between spatial patterns of annual mean ice thickness 24 bias and ice thickness seasonal cycle amplification. However, the external forcings (melt onset and downwelling 25 LW bias) cannot entirely explain the spatial patterns in the mean sea ice state biases, because on a regional scale 26 effects of sea ice convergence, and hence dynamics, become more important. The annual mean ice thickness 27 bias seen in HadGEM2-ES is associated with a thickness maximum on the Pacific side of the Arctic, at variance 28 with observations which show a similar maximum on the Atlantic side. It was shown by Tsamados et al (2013) 29 that such a bias could be reduced by introducing a more realistic sea ice rheology.   The ISF analysis, as presented, does not comprise an exhaustive list of processes affecting Arctic Ocean surface 28 fluxes. The missing process of the effect of snow fraction on surface albedo has already been noted, and its 29 likely effect on the total June and July ISF bias. We note also that the direct effect of thinning ice on ice albedo 30 could induce an additional flux bias relative to the real world, despite the fact that this effect is not represented 31 in HadGEM2-ES. The effect of snow thickness bias on winter conduction and surface temperature is another 32 such process which cannot be included due to inadequate observations. Model biases in the turbulent fluxes may 33 also be significant. While the process which is likely most important in determining these during the winter is 34 captured (ice fraction in the freezing season), a more detailed treatment of turbulent fluxes would also examine 35 the effect on these of the overlying atmospheric conditions. It is also noted that snowfall itself is a component of the surface flux that could, in theory, be evaluated directly given a sufficiently reliable observational reference.

25
A framework has been designed (the ISF framework) that allows the proximate causes of biases in sea ice 26 volume balance to be separated and quantified. Given reference datasets for independent variables, fields of 27 induced surface flux bias can be calculated from the underlying model bias; these in theory sum to the total 28 surface flux bias. In practice, the total ISF bias matches both the net radiation bias, and the ice volume balance 29 bias to first order: processes evaluated cause around 40cm additional ice growth during the ice freezing season, 30 and 20cm additional ice melt in winter; a missing process for which we have no reference (snow thickness) is 31 likely to account for at least some additional ice melt in summer. However, observational uncertainty in the 32 evaluated terms prevents direct evaluation of the total ISF bias, and is the largest contribution to ISF uncertainty.

33
The ISF analysis enables model biases in sea ice growth and melt rate to be attributed in detail to different 34 causes. In particular, the roles played by the sea ice albedo feedback, by the sea ice thickness-growth feedback, to make model ice thickness both low in the annual mean, and too amplified in the seasonal cycle, with the downwelling LW bias acting to mitigate both effects. The result is consistent with the prediction of DeWeaver 1 et al (2008) that sea ice state is more sensitive to surface forcing during the ice melt season than during the ice 2 freeze season. The analysis also suggests that through an indirect effect on surface albedo at a time when sea ice 3 is particularly sensitive to surface radiation biases, the zero-layer approximation, which was until recently 4 commonplace in coupled models, may be of first-order importance in the sea ice state bias of HadGEM2-ES.

5
The ISF analysis also allows more detailed analysis of the spatial patterns in sea ice volume balance simulation.

6
In particular, the mechanisms behind the near-identical spatial pattern of biases in annual mean ice thickness 7 (likely driven by ice dynamics) and that of biases in the ice volume balance are explicitly demonstrated. Where 8 ice thickness is biased low in the annual mean, an enhanced seasonal cycle is apparent. This is due to the ice 9 thickness ISF bias (in freezing season) and the ice area ISF bias (in melting season), corresponding to the thickness-growth and ice albedo feedbacks. The downwelling LW and melt onset biases, by contrast, are more 11 spatially uniform, and do not contribute to the annual mean ice thickness control on the ice volume balance.

12
The finding that observational uncertainty is the most important cause of uncertainty in the ISF bias calculation 13 itself suggests that if observational uncertainty could be reduced, the ISF analysis could become a very powerful 14 tool for Arctic sea ice evaluation. In particular, large observational uncertainties for snow cover and summer

21
The ISF analysis as presented here is designed specifically to approximate HadGEM2-ES, but could in principle 22 be generalised to other models, particularly by altering the surface albedo parameterisation used here, or by 23 using different sea ice thickness categories. The zero-layer thermodynamic assumption used in the ISF analysis 24 is likely to be appropriate for any model during the ice freezing season, as the largest ISF biases tend to arise 25 from the thinnest ice categories, where the zero-layer approximation is closest to reality. However, there is a 26 question as to whether the zero-layer approximation conceals significant surface flux bias relating to ice sensible 27 heat uptake in the late spring.

28
In the case study presented here, the analysis provides mechanisms behind a model bias in sea ice simulation.

29
However, the analysis could also be used to investigate a sea ice simulation that was ostensibly more consistent 30 with observations, to determine whether or not the correct simulation was the consequence of model biases that conductive flux from the ice surface, unless surface melting is taking place; as the freezing case is being 17 discussed, melting is assumed to be zero. We also make the zero-layer approximation used in HadGEM2-ES, 18 that the sea ice has no sensible heat capacity and that conduction is therefore uniform in the vertical for each

24
Following this, we iterate through the categories, identifying grid cells where the bias is such that a negative 25 category sea ice thickness in the reference dataset is implied; in these cells, the bias is reduced such that the

17
Examining first the winter underestimation (demonstrated in Figure B1 a-c), it is found that for each model

23
However, this effect would have no direct impact on the ISF biases because these are computed from monthly 24 means of the model bias in one variable by the model mean in the other; hence, it is covariance between bias 25 and mean that would induce inaccuracy in this case. By similarly approximating the trend in monthly mean 26 model bias as half the difference between model bias in the adjacent months, the error in downwelling SW and 27 ice area contributions were evaluated. Error in the downwelling SW term was found to be significant early in the 28 summer, with an error of -2.7 Wm -2 in June; error in the ice area term was found to be significant later in the 29 summer, with errors of -1.7 Wm -2 and -1.6 Wm -2 in July and August respectively. However, the August error is 30 small relative to the total ISF bias identified.

3
This apparent problem can be resolved partly by viewing the ISF method as a way not simply of estimating 4 model biases due to a particular variable, but of characterising them, i.e. by accepting that the quantity that we 5 are trying to estimate is itself somewhat subjective. Instead of requiring the ISF method to be correct, it is 6 required that it gives useful, physically realistic results. In the case given above, a sufficient condition is that , a term that can be calculated relatively 11 easily as many of the derivatives go to zero. Averaged over the Arctic Ocean this term was small (below 1 Wm -2 12 in magnitude) in most months of the year, but of significant size in October (3.6 Wm -2 ), due to co-location of 13 substantial negative biases in downwelling LW and category 1 ice thickness in this month, indicating that the 14 true surface flux bias in this month may be substantially smaller (in absolute terms) than the -11.5 Wm -2 15 obtained from summing the ISF biases.