Khaldoun Saleh

Yann H Kerr

Philippe Richaume

Maria-José J Escorihuela

R Panciera

S Delwart

Gérard Boulet

Philippe Maisongrande

J P Walker

P Wursteisen

J P Wigneron

Soil moisture retrievals at L-band using a two-step inversion approach (COSMOS/NAFE'05 Experiment)

Keywords:

des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

Introduction

The European Space Agency will launch the SMOS (Soil Moisture and Ocean Salinity) satellite during the summer of 2009. The generation of land products from SMOS brightness temperature (TB) measurements relies on the inversion of the microwave forward model L-MEB (L-band Microwave model of the Biosphere; [START_REF] Wigneron | L-band Microwave Emission of the Biosphere (L-MEB) model: Description and calibration against experimental data sets over crop fields[END_REF][START_REF] Kerr | SMOS level 2 processor soil moisture algorithm theoretical basis document[END_REF]. The L-MEB model assembles a set of equations describing the emission and scattering of the surface, vegetation and atmosphere at L-band (1.4 GHz). The optimisation of L-MEB for different surfaces has been addressed by numerous studies over the last twenty years, most of them based on the analysis of ground-based L-band data. This paper is part of the COSMOS ESA study (Saleh et al., 2007a), designed to test L-MEB over crops and prairies using airborne L-band data from the COSMOS/NAFE campaign [START_REF] Panciera | The NAFE'05/CoSMOS data set: Toward SMOS soil moisture retrieval, downscaling, and assimilation[END_REF]. While the main focus of the COSMOS study is to address open issues regarding the microwave modelling of land surfaces, dedicated campaigns are planned with specific focus on the validation of land products [START_REF] Delwart | SMOS validation and the COSMOS campaigns[END_REF]. A key aspect concerns modelling of the soil emission, as the relationship between the soil emissivity and the surface soil moisture (SM) at L-band is influenced by the dielectric profile at and near the surface as well as by the effective temperature. As a consequence, a sitespecific characterisation of the soil emission, which is usually represented by the roughness parameters H R [START_REF] Choudhury | Effect of surface roughness on the microwave emission from soils[END_REF], and N R,p [START_REF] Wang | Multi-frequency measurements of the effects of soil moisture, soil texture, and surface roughness[END_REF], is often conducted in soil moisture retrieval studies. These parameters are used to express the ratio between the reflectivity of a rough and a smooth surface through the exponential term exp (-H R • cos NR,p (θ)), with N R,p being dependent on the measured polarisation. The main concern regarding roughness is its effect on soil moisture retrievals. If roughness is underestimated so will be the modelled surface emissivity, and the retrieved soil moisture will be underestimated as a result.

While several studies have addressed the generalisation of H R as a function of the measured surface roughness [START_REF] Mo | A parameterization of the effect of surface roughness on microwave emission[END_REF][START_REF] Wigneron | A simple parameterization of the L-band microwave emission from rough agricultural soils[END_REF] a global roughness correction is not available, and alternative methods need to be developed should roughness also play a role at satellite low-resolution scales (40 km). This is of great importance as rather high microwave roughness (H R N 0.3) was measured in experimental plots observed by ground radiometers at L-band [START_REF] Escorihuela | A simple model of the bare soil microwave emission at L-band[END_REF]Saleh et al., 2007b;[START_REF] Wigneron | L-band Microwave Emission of the Biosphere (L-MEB) model: Description and calibration against experimental data sets over crop fields[END_REF]. These values are comparable to estimates from satellite data at higher microwave frequencies [START_REF] Pellarin | Global simulation of brightness temperature at 6.6 and 10.7 GHz over land based on SMMR data set analysis[END_REF], more sensitive to surface roughness. Therefore, L-band roughness could play an important role in soil moisture retrievals from satellite, and currently it is unclear whether its effect will be reduced by the lower spatial resolution and lower radiometric sensitivity of satellite radiometers compared to ground radiometers.

A study on the validation of the L-MEB model from high-resolution airborne data of the COSMOS/NAFE campaign was recently completed by [START_REF] Panciera | Evaluation of the SMOS L-MEB passive microwave soil moisture retrieval algorithm[END_REF]. The study suggested that while the vegetation parameters available in the literature are generally well suited to simulate the vegetation optical depth (τ NAD ), knowledge of soil roughness is essential. The authors estimated H R to be dependent of soil moisture, with H R ~1 for clay soils and soil moisture 0.20 m 3 m -3 . In order to determine soil roughness, the authors conducted site-specific retrievals of H R requiring measurements of the soil moisture and vegetation water content.

In this paper we propose an alternative method to retrieve soil moisture in areas of unknown roughness and unknown vegetation water content with a view on operational applications. The method is based on a two-step two-parameter (2-P) inversion of the L-MEB model. First, H R and the microwave optical depth at nadir are retrieved in small and homogeneous areas of known soil moisture, and a roughness correction dependent on land use only is obtained (either for crops or for prairies). This roughness correction is then applied to 2-P retrievals of soil moisture and optical depth beyond the calibration area, and it is shown to be valid at different footprint resolutions including mixtures of land use across the whole soil moisture range. The two-step inversion method was previously applied to ground-based L-band measurements in a natural grass plot (Saleh et al., 2007b), where the 2-P inversion of τ NAD and H R provided the best approach to retrieve soil moisture in a surface with organic debris. However, the approach had not been tested before using L-band measurements acquired at different resolutions and over different surfaces. To simulate the brightness temperatures of the surface, the study uses the SMOS ground processor breadboards. The SMOS breadboards are a set of routines using the same microwave model and parameterisations included in the SMOS Level 2 processor [START_REF] Kerr | SMOS level 2 processor soil moisture algorithm theoretical basis document[END_REF], as well as the same rationale concerning the aggregation of simulated brightness temperatures over mixed pixels.

This paper starts with a description of the experimental area and data sets, followed by the methodology involved in processing airborne passive microwave measurements. The parameter retrieval method is described afterwards, followed by the discussion of soil moisture retrievals over crops, grass, and mixed land uses at different resolutions.

Data set

The Goulburn River catchment

The area of study is within the Goulburn River catchment (31.77 S to 32.85 S, 149.67 E to 150.60 E) in New South Wales (Australia). This is a semi-arid region of approximately 6540 km 2 , with a central plateau dedicated to crop growth and grazing. The climate of the catchment is sub-humid to temperate, with annual rainfall of approximately 700 mm in 2005 [START_REF] Rudiger | The Goulburn River experimental catchment data set[END_REF]. The experimental area during COSMOS/NAFE covered four focus farms in the Krui river sub-catchment, on the northwest side (562 km 2 ), and four in the Merriwa river sub-catchment, on the northeast side (651 km 2 ). A location map is provided in Fig. 1a, while Table 1 summarises the characteristics of each focus farm used for this study in terms of soil type and vegetation. Radiometer calibration tests during COSMOS/ NAFE took place at the Glenbawn Lake (Fig. 1a), east of the Goulburn River catchment, and also at the Alexandrina Lake near Adelaide (35°25′S 139°07′E) for COSMOS flights.

Airborne measurements

EMIRAD

COSMOS flights used the EMIRAD fully polarimetric radiometer (Rotbøl et al., 2003) onboard an Aero Commander 500S Shrike aircraft. Onboard instrumentation also included an inertial navigation unit combined with a GPS receiver for aircraft attitude and position recording. Two Potter horn antennas were installed along-track (AT), one facing nadir (37.6°half-power beam-width), the other at an angle of 40°towards the rear of the aircraft (30.6°half-power beam-width). Microwave data were produced at 8 ms integration time. Baseline calibration was performed in the laboratory using a Standard Noise Generator (LN2) and an internal hot load, as well as in the aircraft to account for antenna cable losses. Once airborne, internal calibrations were performed to adjust the gain and noise temperature for each measurement. The measured radiometric performance was 0.9 K for 8 ms integration time and 300 K input temperature. The whole EMIRAD data set covered six weeks of data (14 November-9 December 2005), and approximately 30 h of flight time. However, in this study we only made use of EMIRAD flights supported by ground soil moisture measurements at high resolution (3 flights over grass at Stanley, Roscommon and Midlothian, and 3 flights over crops at Pembroke, Illogan and Merriwa Park, Table 2). These flights were single lines performed over each farm with a nominal footprint size of 350 m conducted in the early morning between 15 November and 23 November. Potential instrumental errors in the rear antenna at horizontal polarisation are being investigated and have been discarded for this study. Therefore, the analysis of EMIRAD data is based on nadir acquisitions at horizontal (H) and vertical (V) polarisation, and 40-degree observations at V polarisation. As a result, EMIRAD data were not used to retrieve soil moisture but instead they were used to either compare forward model simulations to brightness temperature measurements, or for the retrieval of optical depth. Finally, infrared measurements at nadir were obtained concurrently to L-band data. For that purpose a Heimann KT15 sensor was used (full beam-width is 4°).

PLMR

The Polarimetric L-band Multibeam Radiometer (PLMR) was installed onboard an Eco-Dimonna aircraft. The radiometer measured horizontally and vertically polarised radiation at six angles of incidence (+/-7°, +/-22°, +/-38.5°) through an antenna patch array with half-power beam-width between 13°and 16.5°. Modifications to the manufacturer calibration were based on external measurements of the sky brightness temperature (cold point) and that of a microwave absorbent (warm point). These measurements were performed at ground level before and after each flight. Further details are available in [START_REF] Panciera | The NAFE'05/CoSMOS data set: Toward SMOS soil moisture retrieval, downscaling, and assimilation[END_REF], where the complete PLMR flight plan can be found. The PLMR data set covered four weeks (30 Oct-25 Nov 2005), with the last two weeks overlapping with EMIRAD flights. The flights considered for this study are also summarised in Table 2. These included six flights in the along-track direction (AT-PLMR) with nominal nadir footprint of 250 m at -3 dB; flights in the across-track direction with a nominal footprint size at nadir equal to 60 m (CT1- PLMR), and 250 m (CT2-PLMR), all supported once a week by ground sampling at each farm. PLMR data were used to perform retrievals at all farms except at Illogan, where the number of multi-angular footprints was considered to be insufficient.

s = 0.2, c = 0.3 1.5-2.5 Roscommon (14) Grazing s = 0.7, c = 0.1 1.5-2 Midlothian (22) Grazing s = 0.1, c = 0.7 1.5-2 Merriwa Park (23) Wheat s = 0.2, c = 0.4 1.5-2 Cullingral (24) Wheat s = 0.2, c = 0.4 1-2

Ground measurements and ancillary data

The ground data used in this study is a subset of the COSMOS/ NAFE'05 data set [START_REF] Panciera | The NAFE'05/CoSMOS data set: Toward SMOS soil moisture retrieval, downscaling, and assimilation[END_REF] and it comprises the following measurements:

-Soil moisture. Surface soil moisture (SM) measurements were obtained from mobile probes (HydraProbe and ThetaProbe). These are capacitance probes with sampling depth of approximately 6 cm. Sampling spatial resolutions varied between 6 m and 500 m, and for this study an area of approximately 1 km 2 within each farm labelled as HRES2 was considered (Fig. 1b). A smaller area of about 150 m × 150 m within HRES2 included 6-m and 12-m sampling only, and we will refer to it as HRES1 (Fig. 1b). All flights analysed in this study were supported by soil moisture samples both in the HRES1 and HRES2 areas. Soil moisture sampling took approximately 5 h, while the duration of flights varied between 15 min and 2 h. However, we found that temporal variations of soil moisture were less significant than spatial ones (not shown), and the whole sampling period was used to calculate ground-truth soil moistures. -Soil temperature. Surface temperature (T sfce ) at 2.5 cm depth was obtained from temperature monitoring stations located within each farm. Deep soil temperature (T depth ) was obtained from one station in the Krui area measuring soil temperature at 60-cm depth. -Soil texture. Soil texture was estimated from 30-cm depth soil samples obtained in the HRES1 area at each focus farm, with the exception of Cullingral and Roscommon where texture was obtained from soil samples within 5 cm of the surface.

In addition, the following ancillary maps were used to characterise the surface land use and topography:

-Land use. A Landsat 5TM-derived land use supervised classification was used to determine land uses covered by each footprint. The classification was for October 2005 and the resolution was 25 m. The classification provided thirteen categories, grouped for this study into four: prairies, forests, open woodland, and crops. -Leaf area index (LAI). Leaf area index maps were derived from MODIS data [START_REF] Duchemin | Monitoring wheat phenology and irrigation in Central Morocco: On the use of relationship between evapotranspiration, crops coefficients, leaf area index and remotely-sensed vegetation indices[END_REF] at a resolution of 250 m on days 02/11, 11/11 and 17/11. -Digital elevation model. A 250-m digital elevation model of the catchment area was also considered for the processing of airborne data described next.

Methodology

Aircraft data processing

The retrieval of surface parameters like soil moisture is based on the comparison between simulations and airborne measurements of the surface brightness temperature. For these two brightness temperatures to be comparable a number of processing steps are required, as summarised next.

Step 1. Ground location of the antenna footprint

For this step, a regular grid with the antenna directional cosines (function of zenith and azimuth angles in the antenna frame) is projected onto the surface in order to calculate the location and local incidence angle of each grid point [START_REF] Tenerelli | Earth-viewing L-band radiometer sensing of sea surface scattered celestial sky radiation-Part I: General characteristics[END_REF]. This step takes into account the orientation of each antenna beam, and the aircraft attitude corresponding to each measurement.

Step 2. TB simulations at the surface level

This step is used to retrieve necessary ground information (e.g. land use) at each point of the projected grid. Simulated brightness temperatures of the whole area covered by each antenna are calculated using the SMOS Level 2 breadboards [START_REF] Kerr | SMOS level 2 processor soil moisture algorithm theoretical basis document[END_REF] at each grid point.

Step 3. TB simulations at the antenna level

For each acquisition, the rotation of the polarisation basis from the surface to the antenna level due to the aircraft movement is calculated, in order to obtain brightness temperatures at the antenna level from brightness temperatures at the surface. For this rotation (Eq. ( 11) in [START_REF] Waldteufel | About off-axis radiometric polarimetric measurements[END_REF], the 3rd and 4th Stokes parameters are assumed to be zero at the surface. Brightness temperatures at the antenna level are then weighted by the antenna gain, and in this way become ready for comparison with the radiometric measurements.

Given the characteristics of the site we considered footprints to be flat. Footprints including altitude variations higher than 5 m across the -3 dB area were discarded. Note that Steps 1 and 3 could be generally neglected for land applications for stable flights (pitch and roll under 2°), small antenna apertures (if modelling assumes that all radiation originates in the boresight direction) and homogeneous surfaces. They were included in the processing though to account for the multi-angular signature integrated by large aperture antennas, and for a better radiometric description of heterogeneous areas.

L-MEB main equations

The L-MEB model main equations for soil and vegetation are summarised below. Symbols p and θ indicate a dependence on polarisation and angle respectively. The soil brightness temperature TB S ( 1 et. al., 2007).

TB S;p;θ = 1-r 4 S;p;θ exp -H R cos θ ð Þ N R;p T S = 1-r S;p;θ T S ð1Þ
T S = T depth + T sfce -T depth SM w 0 b w 0 : ð2Þ
The final equation for the surface brightness temperature TB (3) accounts also for vegetation, and a soil-vegetation composite temperature T GC described in [START_REF] Wigneron | L-band Microwave Emission of the Biosphere (L-MEB) model: Description and calibration against experimental data sets over crop fields[END_REF]. In this study we considered vegetation and surface temperatures to be equal based on the data availability. Finally, the vegetation parameters in the model are the single scattering albedo ω p , the transmissivity γ and the optical depth parameters τ NAD and tt p .

TB p;θ = 1-ω p 1 -γ p;θ 1 + γ p;θ r S;p;θ + 1-r S;p;θ γ p;θ T GC ð3Þ
γ p;θ = exp -τ NAD sin 2 θ ð Þ + tt p cos 2 θ ð Þ = cos θ ð Þ: ð4Þ

Parameter retrievals

The retrieval of model parameters and surface properties is based on the optimisation of the L-MEB model when contrasted with measurements as described in Section 3.1. The inversion technique uses a modified least squares cost function (5) with the possibility to Author's personal copy constrain the retrieved parameters (p i ) and brightness temperatures (measured TB°, and modelled TB).

C F = P TB B θ;p -TB θ;p 2 σ 2 TB + X i p i -p ini i 2 σ 2 p : ð5Þ
In this study retrieved parameters H R , SM, and τ NAD were allowed to change with variance σ p 2 = 1. Such variance is higher than that expected one at the field scale as shown by field measurements, and therefore represents a low constraint on the retrieved parameters. The standard deviation for TB measurements σ TB was equal to 1 K. Either H R and τ NAD , or SM and τ NAD was retrieved from multiangular dual-polarisation PLMR data; for EMIRAD, only two configurations were available (nadir, and 40°at vertical polarisation), and either direct simulations of the brightness temperature or retrievals of one parameter (τ NAD ) were conducted. In the latter case, one optical depth was retrieved from simulations of TB at nadir and at 40°at vertical polarisation. The initial values for the retrieved parameters were set to 0.1 m 3 m -3 for SM, 0.1 for τ NAD and 0.1 for H R .

Finally, for the selection of multi-angular TB measurements, the surface was gridded at resolutions between 250-m and 500-m as detailed later in the paper. Footprints with the beam centre included in each cell were used to retrieve parameters that were taken as representative of each cell. Standard deviations of each retrieved parameter were also calculated as part of the optimisation procedure.

Results

In the soil moisture sampling areas both radiometers values were well above those simulated for a surface with small roughness, moderate biomass, and the measured soil moisture (e.g. H R =0.1, τ NAD = 0.15). Both the increased level of TB with respect to a planar surface as well as the low angular ratio between nadir and off-nadir measurements could only be explained in terms of the surface roughness through the roughness parameters H R and N R,p (1). This feature could not be attributed to instrumental errors in the aircraft data, ground instruments, surface temperature, or soil moisture, and was observed at different flying altitudes and both radiometers. The approach presented next aimed at minimising the effect of roughness on soil moisture retrievals using microwave measurements and knowledge of the land use only.

Soil emission characterisation

Estimates of soil parameters were obtained from CT1-PLMR flights (Table 2) in the homogeneous HRES1 areas, where a large number of soil moisture measurements were available (typically above 200). The reference soil moisture for the HRES1 area was the simple average of all soil moisture measurements within that area (SM field ). As summarised in the Introduction, the selected approach was a two-parameter (2-P) inversion that fixes soil moisture, and estimates simultaneously H R and the optical depth at nadir (τ NAD ) from multi-angular brightness temperatures. In fact, 3-parameter (3-P) inversion tests where SM, H R and τ NAD,3P were obtained simultaneously from CT1-PLMR data failed at retrieving SM correctly but provided optical depths comparable to those obtained from the 2-P approach (τ NAD,2P ). While the root mean square error (RMSE) between τ NAD,2P and τ NAD,3P was under 0.04 for grass and crops, the RMSE between SM field and SM 3P was above 0.2 m 3 m -3 whereas errors under 0.05 m 3 m -3 are desirable for satellite products [START_REF] Walker | Requirements of a global near-surface soil moisture satellite mission: Accuracy, repeat time, and spatial resolution[END_REF]. These results suggest that 3-P retrievals fail because the inversion returns low soil moisture in rough soils, as for most configurations these surfaces are characterised by high brightness temperature and low angular ratio at a given polarisation. Therefore, roughness needs to be highly constrained in order to obtain good retrievals of soil moisture as found by other studies [START_REF] Pardé | N-parameter retrievals from L-band microwave observations acquired over a variety of crop fields[END_REF][START_REF] Wigneron | L-band Microwave Emission of the Biosphere (L-MEB) model: Description and calibration against experimental data sets over crop fields[END_REF]. The 2-P inversion with known SM in the HRES1 area produced very good fits in terms of brightness temperatures (RMSE b 3 K, R 2 N 0.95), and clear differences between H R over crops and grass (Fig. 2). Other parameters were tt p =1, ω p =0, N R,H =1 and N R,V =0, which resulted from best-fit model simulations conducted farm to farm with ancillary vegetation information (not discussed here). Based on Fig. 2, average values of H R over crops and grass were obtained.

H R crop ð Þ= 1:0 σ = 0:2 ð Þ H R grass ð Þ= 0:4 σ = 0:2 ð Þ : ð6Þ
Note that at Roscommon farm ('14' in Fig. 2) large roughness was obtained for rather dry soil (0.04 m 3 m -3 ) as well as very low optical Fig. 2. Surface roughness parameter H R retrieved over crops (inside ellipse) and grass (outside ellipse) from CT1-PLMR flights in the HRES1 area. Numbers correspond to farms (labels in Table 1) and bars are standard deviations of averaged field soil moisture in the HRES1 area (SM field , x axis) and retrieved roughness (H R , y axis) respectively. depth compared to previous days. Because field and satellite measurements indicated no significant vegetation changes, we did not include that point in the averaged H R (6). We also determined that differences in the averaged roughness parameter H R were negligible between using a small standard deviation on the optical depth (σ τ = 0.01 in Eq. ( 5)), or a large one (σ τ = 1). The error in soil moisture retrievals in the roughness calibration areas (HRES1) was RMSE = 0.040 m 3 m -3 (R 2 = 0.94) for grass, and RMSE = 0.046 (R 2 = 0.94) m 3 m -3 for crops (Table 3), following a 2-P inversion of SM and τ NAD and the constant roughness values in Eq. ( 6).

The same roughness fit was applied to PLMR along-track data over crops (AT-PLMR, Table 2) in the HRES1 areas. These flights were obtained at a higher altitude (250 m nominal resolution), and the accuracy of the soil moisture retrievals was RMSE = 0.044 m 3 m -3 (R 2 = 0.68).

Alternative methods to determine surface roughness were also explored, such as "field-to-field" 1-parameter (1-P) retrievals of H R with known soil moisture, and optical depth estimates based on field measurements of vegetation water content and satellite-derived LAI. Results of the 1-P inversion of H R led to mean roughness H R = 1.1 for crops, and mean H R = 0.5 for grass for τ NAD =0.06•LAI, with 0.06 being an optimised parameter. These roughness values are close to those derived from the 2-P approach in Eq. ( 6), the main difference between the two being a higher sensitivity of H R to SM observed in the 1-P estimates of H R for crops. For example, in the 1-P case the roughness variation between dry (0.05 m 3 m -3 ) and wet (0.40 m 3 m -3 ) soils was approximately ΔH R = 0.4 for crops (mean H R = 1.1), and the correlation between H R and SM was R 2 = 0.54. However, applying a constant roughness correction was explored further as a way to retrieve soil moisture in areas of unknown roughness and unknown optical depth as described next.

Soil moisture retrievals in the HRES2 areas

Soil moisture retrievals beyond the calibration area were investigated from low altitude flights (CT1-PLMR), and from higher altitude flights (CT2-PLMR, AT-PLMR and AT-EMIRAD). Reference soil moisture for validation in this region (HRES2) was obtained at either 250-m cells (CT1-PLMR, AT-PLMR, AT-EMIRAD flights) or 500-m cells (CT2-PLMR). The size of the grids was selected according to the size of the footprints, and to ensure multi-angular measurements inside each cell. Within each of these cells, one reference soil moisture was obtained from averaging field soil moisture data, and only cells with more than 20 soil moisture sampling points were compared to retrievals. Therefore, each cell had a corresponding retrieved soil moisture (SM R ) obtained from multi-angular TB measurements centred in the cell, and obtained concurrently to the optical depth τ NAD (Figs. 3 and5). Table 3 summarises the inversion results following this approach, i.e. 2-P inversions of SM and τ NAD with the L-MEB parameters summarised in Table 4. Table 3 also shows results of a 3-P inversion (SM, τ NAD , H R ) for comparison. Results for crops, grass and mixed areas are discussed next in detail.

Crops

Soil moisture retrievals over crops were conducted in 20 cells (including several days) where the crop fraction covered by the footprints was higher than 95%, in order to analyse 'pure' pixels (Fig. 3a). Using the roughness fits in Eq. ( 6) the accuracy in the soil moisture retrievals from CT1-PLMR data was high (RMSE =0.046 m 3 m -3 ), showing that the roughness relation derived from high-resolution flights was still valid outside the calibration area.

As for flights in the HRES1 area, poor correlation (R 2 = 0.3) was found between τ NAD and estimates of LAI based on three MODIS images, the average bias being equivalent to LAI =1 m 2 m -2 for τ NAD =0.06•LAI. However, we found that this is not a limitation for the accuracy of the soil moisture product. For comparison, note that results of 1-P inversions of soil moisture with optimised optical depth (τ NAD =0.06•LAI) had a RMSE = 0.041 m 3 m -3 . 1-P inversions used a roughness fit obtained from individual H R retrievals in areas with known soil moisture (H R =1.4-1.1 • SM, R 2 = 0.54).

PLMR-based results were then compared to AT-EMIRAD flights (TB 0,H , TB 0,V , and TB 40,V ) through direct modelling of the brightness temperature in the HRES2 area (Fig. 4a). Direct simulations for these three configurations covered three different days and included three crop farms (SM range from 0.13 m 3 m -3 to 0.29 m 3 m -3 ). For these simulations a value of τ NAD = 0.13 was used as derived from the PLMR data set (Table 3), mean soil moisture values measured in the HRES2 area, and the roughness fit in Eq. ( 6). For EMIRAD flights we used the surface infrared temperature onboard the aircraft as a proxy for T surf in Eq. ( 2). Thermal infrared measurements had a standard deviation lower than 1 K across each farm, and therefore a constant value was used for the nadir and off-nadir simulations. Based on this approach, the RMSE between simulated and measured AT-EMIRAD brightness temperatures was 2.2 K (R 2 = 0.98, 157 footprints), including simulations at the Illogan farm that were not considered previously for the roughness calibration from PLMR data.

Grass

Soil moisture retrievals from CT1-PLMR flights over grass were possible in 26 cells across three different farms. The same 2-P approach as for crops was followed, and soil moisture could be estimated with an accuracy of 0.034 m 3 m -3 (Fig. 3b). In addition, CT2-PLMR flights over 'pure grass' were analysed. At this resolution (250 m nadir) estimates were still good (RMSE = 0.042 m 3 m -3 , Fig. 5) despite the increased cell size (500-m, in order to obtain multiangular footprints).

In terms of optical depths obtained from CT1-PLMR flights we observed differences between the farms, with mean τ NAD ~0.28 at Stanley and Midlothian, and mean τ NAD ~0.11 at Roscommon where the vegetation biomass was significantly shorter than at other farms in the HRES2 area.

These optical depth values were used for comparison with AT-EMIRAD flights over the same farms (Fig. 4b). In this exercise τ NAD was retrieved for comparison with previous estimates obtained from PLMR data. Very good agreement was found between measured and simulated TBs at the Stanley farm (RMSE = 1.7 K, R 2 = 0.98), where the retrieved τ NAD was 0.30. At Roscommon, the retrieved τ NAD was 0.08 (RMSE = 1.8 K), which is consistent with the low values obtained from PLMR data at this farm too. This agreement in optical depths for fixed roughness and soil moisture indicates that measurements from the two sensors were similar. Results at Midlothian were less conclusive, as soil moisture measurements in the densely sampled areas were highly variable (σ(SM) / mean(SM) ~1) on the day of EMIRAD flights. Outside the dense sampling region retrievals of τ NAD varied between 0.2 and 0.3 (RMSE b 1 K), but given the smaller number of soil moisture points these have not been included in Fig. 3b.

Based on positive results of five of the six farms with AT-EMIRAD flights we consider the data from the two sensors to be in reasonable agreement. Also, the fact that roughness correction appears appropriate at different times of the day (early in the morning for EMIRAD, and up to the early afternoon for PLMR) indicates that surface temperature does not have a significant impact on the roughness correction.

Land use mixing

Finally, footprints including different land uses were analysed. For this exercise PLMR along-track flights were most appropriate. The retrieved parameters were obtained for 250-m cells, and the mixing factor was taken as a cell-averaged value of the land use fractions of individual footprints falling into each cell. Cells at Midlothian and Merriwa Park satisfied the requirements to include multi-angular measurements and multiple land uses (crop-grass mixtures, or grassopen woodland mixtures). An illustration of these cells is shown in Fig. 6 together with an aerial photograph of a woodland area at Midlothian. At Merriwa Park, three cells at the boundary of crop and grass fields were examined across four different days. Note that not all flights provided TB measurements of all cells, although each cell could be monitored at least twice. In terms of mixing, the crop fraction varied between approximately 20% and 60% in the antenna -3 dB region (Table 5). At Midlothian, only areas with a small fraction of open woodland were supported by enough soil moisture sampling, and four cells were selected for the study. The woodland fraction varied between 7% and 30% (Table 5).

Fig. 1 .

 1 Fig. 1. (a) The Goulburn River catchment experimental site: catchments, focus farms, and calibration lake, (b) Example of soil moisture sampling points at Merriwa Park (black dots), high-resolution soil moisture sampling area used for model calibration (HRES1), and extended area HRES2 used for model validation.

  ) is obtained from the Fresnel coherent reflectivity r S,p * and a set of roughness parameters H R and N R,p for non-smooth surfaces. The soil temperature T S (2) considered to contribute to the soil emission is a combination of temperatures, deep and near the surface, modified by a factor dependent on the surface soil moisture (SM), and two fixed parameters w 0 and b w 0 (default values are w 0 =0.3 m 3 m -3 and b w 0 =0.3, Wigneron

Fig. 3 .

 3 Fig. 3. (a) Soil moisture retrievals (SM R ) over 'pure' crops, and 'pure' grass (b) from CT1-PLMR flights in the HRES2 area based on a 2-P inversion (SM R , τ NAD ) with the roughness fit in Eq. (6). Bars indicate the standard deviation of the field soil moisture (SM field ) and the retrieved soil moisture (SM R ).

Fig. 4 .

 4 Fig. 4. Comparison between AT-EMIRAD measurements (TB 0,H , TB 0,V , TB 40,V ) and simulations of the brightness temperature for crops (a), and grass (b) based on PLMRderived roughness (6). Line 1:1 (solid), and +/-3 K (dashed).

Fig. 5 .

 5 Fig. 5. Soil moisture retrievals from CT2-PLMR flights over grass in the HRES2 area based on a 2-P inversion (SM R , τ NAD ) with the roughness fit in Eq. (6). Bars indicate the standard deviation of the measured field soil moisture (SM field ) and the retrieved soil moisture (SM R ).

Table 1

 1 Focus farms at the Goulburn river catchment.

	Farm	Land use	Sand (s),	LAI-range (m 2 m -2 )
			clay(c) fraction	from MODIS 250 m
	Pembroke (11)	Wheat	s = 0.06, c = 0.7	2-2.5
	Stanley (12)	Grazing	s = 0.06, c = 0.6	1.5-2
	Illogan (13)	Barley, oats

Table 2

 2 Summary of flights supported by intensive ground sampling used for this study.

	Label	Sensor	Configuration Resolution	Land	Number	Number
				-3 dB nadir	use	of farms	of flights
	AT-EMIRAD EMIRAD Along-track	350 m	Crop	3	3
	AT-EMIRAD EMIRAD Along-track	350 m	Grass 3	3
	AT-PLMR	PLMR	Along-track	250 m	Crop	2	4
	CT1-PLMR	PLMR	Across-track	60 m	Crop	3	10
	CT1-PLMR	PLMR	Across-track	60 m	Grass 3	10
	CT2-PLMR	PLMR	Across-track	250 m	Grass 3	9

Table 3

 3 Results of parameter retrievals.

	Flight	Area	Roughness	RMS SM	R 2 SM No. of	τ NAD
				[m 3 m -3 ]		cells
	Soil moisture retrievals -crops
	CT1-PLMR HRES1	3-P inversion N0.1	b0	1 2	0 . 1 4
			2-P, H R = 1.0	0.046	0.94	12	0.13
	AT-PLMR	HRES1	3-P inversion N0.1	0.50	5	0.13
			2-P, H R = 1.0	0.044	0.68	5	0.13
	CT1-PLMR HRES2	3-P inversion N0.1	b0	2 0	0 . 1 2
			2-P, H R = 1.0	0.046	0.93	20	0.13
	Soil moisture retrievals -grass
	CT1-PLMR HRES1	3-P inversion N0.1	b0	8	0 . 1 5
			2-P, H R = 0.4 0.040	0.94	8	0.15
	CT1-PLMR HRES2	3-P inversion N0.1	b0	26	0.25
			2-P, H R = 0.4 0.034	0.96	26	0.25
	CT2-PLMR HRES2	3-P inversion N0.1	b0	18	0.24
			2-P, H R = 0.4 0.042	0.95	18	0.25
	Soil moisture retrievals -mixed ( a H R of dominant land use (LU), b H R function of LU)
	AT-PLMR	HRES2-crop 2-P, a H R	0.059	0.94	9	-
			2-P, b H R	0.046	0.98	9	-
		HRES2-open	2-P, H R = 0.4 0.026	0.70	4	-
		woodland
	TB direct simulations -EMIRAD. (⁎see text)
	Flight	Area	Roughness	RMS TB [K] R 2 TB No. of pts τ NAD ⁎
	AT-EMI	HRES2-crop 2-P, H R = 1.0	2.2	0.98	157	0.13
	AT-EMI	HRES2-grass 2-P, H R = 0.4 1.7	0.97	87	0.08, 0.3

Table 4 L

 4 -MEB parameters used for retrievals.

	K. Saleh et al. / Remote Sensing of Environment 113 (2009) 1304-1312					1309
	Land use	H R	N R,H	N R,V	ω H	ω V	tt p
	Crop	1.0	1	0	0	0	1
	Grass	0.4	1	0	0	0	1
	Open woodland	0.4	1	0	0	0.09	1

K.Saleh et al. / Remote Sensing of Environment 113 (2009) 1304-1312

Acknowledgments

The authors would like to acknowledge the COSMOS and NAFE'05 teams as well as J. Tenerelli and N. Reul (Boost/Ifremer) for sharing their processing tool TRAP (The Terrestrial Radiometry Analysis Package, version 0.1). This study was conducted under the European Space Agency contract RFQ-3-11801/06/NL/FF and the CNES/TOSCA program.

Author's personal copy

The same approach based on a 2-P inversion of soil moisture and τ NAD was followed, and the roughness distribution within the footprint was considered for the simulations. More precisely, the roughness values corresponding to each antenna grid point were chosen based on the land use, this is H R = 1 for crops, and H R = 0.4 for grass and open woodland.

Retrievals in cells including crops and grass at Merriwa Park (Fig. 6a) showed that the mean RMSE between the field soil moisture and the retrieved one was 0.046 m 3 m -3 (R 2 = 0.98) when the roughness of each land use included in the footprint was considered to compute the brightness temperature (Table 5). The error was 0.059 m 3 m -3 (R 2 = 0.94) when the roughness of the dominant land use (crop or prairie) was used instead.

For woodland patches in grass fields (Fig. 6b,c) the same roughness correction as that of grass was assumed, and vertical scattering albedo ω V was 0.09 ('default SMOS Level 2'). Soil moisture retrievals in these cells were close to the field measurements with RMSE = 0.026 m 3 m -3 (R 2 = 0.70). Because the fraction of forest vegetation in open woodland areas was small, using the same roughness correction for grass and open woodland agreed with the observations. These results are positive as they provide good prospects for retrievals in mixed land uses, and highlight the importance of including appropriate roughness corrections in the simulation of mixed footprints.

Summary and discussion

This paper investigated L-band measurements over grass and crops obtained from two airborne-based L-band radiometers (PLMR and EMIRAD) during the COSMOS/NAFE'05 campaign in south-east Australia. The main objective of the study was to evaluate the performance of the microwave radiative transfer model which sits at the core of ESA's SMOS soil moisture algorithm, and to point out differences between modelling outputs and measurements. In this respect, the analysis of the COSMOS/NAFE'05 data set showed that calibration of the surface roughness is key for successful soil moisture retrievals at L-band, at least at the investigated spatial resolutions (60-375 m nadir), and an approach to overcome the limitations of unknown roughness was presented. Many studies report site-specific calibrations of surface roughness, but we still lack a clear understanding of the physical meaning of such roughness, and how it could be applied globally to minimise its effect on soil moisture retrievals. The approach developed uses a two-step two-parameter inversion of the radiative transfer model (L-MEB). This approach was tested from ground-based L-band measurements at high resolution before (Saleh et al., 2007b), but not at larger scales as observed from an aircraft, and over different types of surfaces.

First, the L-MEB main roughness parameter H R and the optical depth τ NAD were obtained simultaneously in small areas of known soil moisture using multi-angular radiometric brightness temperatures at H and V polarisations (PLMR), and a roughness correction based on land use only was obtained. Then, simultaneous soil moisture SM and optical depth τ NAD retrievals with a high constraint on roughness (effectively a 2-P inversion) were conducted beyond the calibration area, at different resolutions (60 m-250 m), and viewing configurations (across-track, along-track). The soil roughness correction was also used to simulate measurements of the second radiometer (EMIRAD) as a mean of validation, with errors in the simulated brightness temperatures around 2 K. Globally, soil moisture estimates were obtained with errors ranging between 0.03 m 3 m -3 and 0.05 m 3 m -3 compared to the mean field soil moisture, and included 'pure crop' and 'pure grass' areas, as well as mixed areas of crop and grass, and grass and open woodland.

The main emphasis of the two-step two-parameter inversion approach is in reducing the influence of surface roughness to retrieve soil moisture from L-band data. More studies though are needed to evaluate the options for obtaining vegetation products from this approach. The fact that a low correlation between retrieved optical depths and leaf area index was found could be linked to using a constant roughness correction, as soil roughness has been shown to decrease with soil moisture in some studies. Soil moisture has an impact on the wave penetration depth, and roughness estimates could be affected by the thickness of the sensing layer, though no conclusions on this matter have been found yet. Nevertheless, the main interest of the method lies in the opportunity for the calibration of surface roughness should it appear crucial at satellite scales, as only land use knowledge would be required.

From an operational perspective, two steps towards calibrating roughness could be envisaged, first at the instrumented mission's validation sites distributed across the world, then elsewhere. At the validation sites, two-parameter retrievals (H R and τ NAD ) would provide the first estimates of roughness at the scale of a satellite footprint for the whole soil moisture range. In parallel, time series of brightness temperatures could be used to detect very dry and very wet soil to Author's personal copy analyse roughness estimates obtained for extreme soil moisture conditions only. The interest of such experiment would be in investigating the extension of the calibration approach to noninstrumented areas, where very dry or very wet soil conditions could be determined either from brightness temperature time series, or change-detection techniques applied to microwave active measurements [START_REF] Wagner | Method for estimating soil moisture from ERS scatterometer and soil data[END_REF]. In this way, studies pursuing a global correction of roughness based on the land use could be explored. Finally, the continuation of field experiments to improve our understanding of the physical basis for the emission of rough soils at L-band is essential, and upcoming airborne experiments will contribute to better understand the role of scaling in surface roughness.


Sébastien Destercke
email: sebastien.destercke@hds.utc.fr

Agnes Rico
email: agnes.rico@univ-lyon1.fr

Olivier Strauss
email: olivier.strauss@umontpellier.fr

Using atomic bounds to get sub-modular approximations

Set functions and inferences over them play an important role in many AI-related problems. They are commonly used in kernel-based methods, game theory, uncertainty representation and preference modelling. Of particular importance are additive set-functions and their associated expectation operator, as well as sub-modular functions. However, specifying precisely and completely such set functions can be a daunting task as well as a lot of data and knowledge. When time or information is missing, it may be handy to have approximating tools of such set-functions, that require less information and possibly enjoy nice mathematical properties. In this paper, we show that if the only information we have are atomic bounds of such functions, we can build conservative approximations that are either sub-or super-modular. We then illustrate the potential use of our approximation on a convolution-based signal processing problem.

Introduction

Set functions or measures play an important role in many aspects of artificial intelligence. The most well-known are undoubtedly probability measures, i.e., additive, normalised positive set functions. Increasing set functions, called capacities or fuzzy measures, generalise such probability measures. They appeared in 50's, and have been extensively used in areas related to decision [START_REF] Cohen | Cardinal extensions of the eu model based on the choquet integral[END_REF]. In particular, they can represent the weights of coalitions in cooperative games, or measures of uncertainty such as lower probabilities or belief functions (see [START_REF] Grabisch | Set functions, Games and Capacities in Decision-Making[END_REF]). Particular fuzzy measures that are maxitive rather than additive are possibility measures [START_REF] Zadeh | Fuzzy sets as a basis for a theory of possibility[END_REF], the basic building blocks of possibility theory [START_REF] Prade | Possibility Theory. An Approach to Computerized Processing of Uncertainty[END_REF][START_REF] Prade | Possibility theory and its applications: Where do we stand?[END_REF].

Other capacities that have attracted a lot of attention are sub-modular functions, as they have taken a very important role in many aspects of artificial intelligence. They are used in numerous problems of machine learning and signal processing [START_REF] Tohidi | Submodularity in action: From machine learning to signal processing applications[END_REF][START_REF] Victor Gabillon | Adaptive submodular maximization in bandit setting[END_REF], in social choice theory where they can be used to induce fairness in the obtained solutions [START_REF] Benabbou | Interactive optimization of submodular functions under matroid constraints[END_REF], and in modelling preferences in item selection problems [START_REF] Tschiatschek | Selecting sequences of items via submodular maximization[END_REF]. Indeed, the sub-modularity property comes with many advantages, notably the one to be able to use greedy optimisation algorithms with good guaranteed approximation properties.

There are good reasons to consider set functions that are not necessarily increasing with respect to set inclusion, nor that are normalized as probability or possibility measures can be. A question that follows is whether approximation tools allowing for efficient inferences are still valid in such contexts. This is the question we address in this paper, where we focus on approximation obtained by considering atomic bounds (i.e., lower/upper bounds over singletons).

The paper is structured as follows. Section 2 is a brief recall of set-functions and their use. We will focus on additive set functions and sub-modular set functions. Section 3 is devoted to the approximation of set of additive set-functions and sub-modular functions. Section 4 presents an application to signal convolution.

Set-functions and their use

In this paper, we are concerned with approximating (sets) of discrete set-functions. By set function κ over a finite space X = {x 1 , . . . , x n }, we understand a bounded realvalued function κ : 2 X → R from the power set of X to the reals. We furthermore assume that κ(∅) = 0. The conjugate κ c of a set function κ is κ c (A) = κ(X ) -κ(A c ).

Additive set-functions and expectation

One of the most widespread set functions are additive ones, i.e., set functions defined from a mass distribution ψ : X → R assigning a mass to each element, from which we have κ ψ (A) = x∈A ψ(x) (in the rest of the paper, we denote by κ ψ additive measures derived from the mass distribution ψ).

Once we are given such a set-function, a very common operator applied to them to make inferences is the linear expectation operator. If we have a real valued bounded function f : X → R, such an operator is

E κ ψ (f ) = x∈X f (x)ψ(x) (1)
= n i=1 (f (x [i] ) -f (x [i-1] ))κ ψ (A [i] ), (2)
where [.] is the permutation such that (f 2) will be helpful when considering sub-modular functions later on. If we assume that ψ is non-negative, ψ ≥ 0, and that κ ψ (X ) = x∈X ψ(x) = 1, we retrieve classical probabilities and (1) is just standard expectations. Also, ψ being positive ensures that the obtained function κ is monotonic with respect to inclusion, i.e., κ ψ (A) ≥ κ ψ (B) whenever B ⊆ A.

(x [1] ) ≤ • • • ≤ f (x [N ] )), A [i] = {x [i] , . . . , x [N ] } with, by convention, x [0] = 0 and A [N +1] = ∅. Equations (
In this paper, however, we make neither of these two assumptions, i.e., κ ψ (X ) can be any real value, as well as ψ(x) that can be negative. This latter relaxation implies that we may have B ⊆ A but κ ψ (A) > κ ψ (B).

While the latter property may seem odd to some readers, there are good reasons to consider such a relaxation. For instance, deriving a signal through convolution (the application of Section 4 is of this kind), we are considering a non-null function ψ such that κ(X ) = 0, meaning that some values of ψ must be negative. Similarly, in social ranking problems [START_REF] Moretti | Some axiomatic and algorithmic perspectives on the social ranking problem[END_REF], the elements x i can be individuals and one search to assess/order the performances of groups or subsets of individuals. In this case, it can be that two groups A ⊂ B are such that B performs worse than A (A B), for instance because B contains individuals that do not get along together and worsen the whole group performance. Should such a statement be numerically encoded as a set function, we should have κ ψ (A) > κ ψ (B).

Example 1. Consider a digital signal f taking values on two successive points:

f (x 1 ) = 3, f (x 2 ) = 5 computing its derivative in x 2 comes down to compute E κ ψ with ψ(x 1 ) = -1, ψ(x 2 ) = 1, giving E κ ψ =2

Sub-modular set functions

Another kind of set functions that are widely used in artificial intelligence are submodular set functions. Sub-modularity means that, for any pair A, B ⊆ X of event, κ is such that

κ(A ∪ B) + κ(A ∩ B) ≤ κ(A) + κ(B), (3)
or equivalently [15, Ch. 44] that, for any A ⊆ B and any x ∈ B, we have

κ(A ∪ {x}) -κ(A) ≥ κ(B ∪ {x}) -κ(B). (4)
This latter formulation is sometimes referred as the property of diminishing return, modelling the fact that element {x} has a higher value when added to a smaller subset. However, in contrast with the additive case, specifying a sub-modular set-function requires in general to specify or to work with 2 |X | values, a daunting quantity. In practice, it may be useful to have approximations requiring to collect and work with a limited number of values. In this paper, we will see that this can be done by simply considering 2|X | values of κ. Note that the idea of non-monotonic (w.r.t. inclusion), nonadditive set functions already appears in [START_REF] Peter | Binary interactions and subset choice[END_REF], precisely to model numerically rankings of subsets.

Sub-modular set functions have many interesting properties, but we will only mention those useful for our purpose, and proved in [START_REF] De Waegenaere | Choquet integrals with respect to non-monotonic set functions[END_REF]. A first property is that its core

C κ = {ψ : ∀A ⊆ X , κ ψ (A) ≤ κ(A), ψ(X ) = κ(X )}, (5)
that is the set of additive measures that is dominated by it is never empty. In particular, this means that sub-modular set functions can be useful to describe or approximate sets of additive set functions, e.g., in convolution problems where the convolution filter is ill-known.

A second property one can also show [START_REF] De Waegenaere | Choquet integrals with respect to non-monotonic set functions[END_REF] that applying the Choquet integral to a bounded function f : X → R gives

E κ (f ) = n i=1 (f (x [i] ) -f (x [i-1] ))κ(A [i] ), (6)
= sup ψ∈Cκ E κ ψ (f ) (7)
Similar properties apply to the super-modular

κ c E κ (f ) = n i=1 (f (x [i] ) -f (x [i-1] ))κ c (A [i] ), (8)
= inf ψ∈Cκ E κ ψ (f ) (9)
3 Approximating set-functions with atomic bounds

In this section, we provide our main results concerning non-additive set functions, that show that one can build sub-modular approximations by only focusing on atomic bounds of some sets of additive measures or of some initial sub-modular function. This provides a very easy way to get approximations by focusing on only 2n values.

Approximating sets of additive set-functions

We first consider the case where we have a set M of possible additive measures. In signal processing, those may come from the fact that the convoluting kernel has an ill-known bandwidth, while in preference modelling they may represent the doubt an agent has about the strength or rankings of coalitions. M may also be a neighbourhood around some initial distribution ψ 0 , used for instance to make a sensitivity/robustness study around ψ 0 . However, working with general sets M is highly inconvenient, whether they are discrete or continuous. We will now show that we can easily obtain a sub-modular approximation of M. To do this, we will denote by

ψ(x) = inf ψ∈M ψ(x), ψ(x) = sup ψ∈M ψ(x).
Let us now define the two set functions

µ M (A) = max x∈A ψ(x), κ ψ (X ) - x∈A c ψ(x) , (10)
µ M (A) = min x∈A ψ(x), κ ψ (X ) - x∈A c ψ(x) , (11)
with κ ψ (X ) the value for any ψ ∈ M. The subscript µ M denotes that µ originates from a set of additive measures.

Example 2. Let us assume that some group of persons want to pack some food before going on a hike. They can choose between (a)pples, (b)ananas, (c)horizo and (d)anish pastries, so that X = {a, b, c, d}. In order to know what to take (it is impossible to take all four types of food), members of the group are allowed to distribute 10 points between all food types to express their individual preferences. They can also give negative scores (up to -5) if they really dislike a food. We assume that food values are additive. So we know that κ ψ (X ) = 10, and the span of scores given by the group to each type of food is summarised in Table 1. Everyone likes apples as even the lower score ψ is positive (some more than others), but some dislikes bananas (e.g., because of their taste) while some dislikes chorizo (e.g., as they are vegetarian). We can now try to assess the possible values of various packages to compare them. For instance, evaluating the set {a, b, c}, we get

Food (X ) a b c d ψ 5 7 4 4 ψ 1 -3 -3 -2
µ({a, b, c}) = min(5 + 7 + 4, 10 -(-2)) = 12 µ({a, b, c}) = max(1 -3 -3, 10 -4) = 6
so that this package is clearly interesting, for instance compared to {b, c, d} for which µ({b, c, d}) = 9.

We will now study the properties of the set functions given by Equations 10-11. Let us first show that the two set functions are conjugate ones. Proposition 1. For any M, µ M and µ M are conjugate.

Proof.

µ M (A) = µ M (X ) -µ M (A c ) = κ ψ (X ) -max{ x ∈A ψ(x), κ ψ (X ) - x∈A ψ(x)} = min{κ ψ (X ) - x ∈A ψ(x), x∈A ψ(x)} = µ M (A).
We can then show that µ M can be used as a conservative approximation of µ M , by the following proposition.

Proposition 2. We have M ⊆ C µ M .
Proof. The proof simply consists in showing that for any ψ ∈ M, we have ψ(A) ≤ µ(A) for any A ⊆ X . To see this, observe that for any ψ, we have ψ ≥ ψ. This means that x∈A ψ(x) ≥ x∈A ψ(x). Similarly, since ψ ≤ ψ, we have

x ∈A ψ(x) ≤ x ∈A ψ(x), hence κ ψ (X ) -x∈A c ψ(x) ≥ κ ψ (X ) -x∈A c ψ(x) = x∈A ψ(x)
where the last equality comes from ψ being additive and x∈X ψ(x) = κ ψ (X ).

We now demonstrate that µ M is a super-modular set function. Among other things, this ensures that applying Equation [START_REF] Devroye | No empirical probability measure can converge in the total variation sense for all distributions[END_REF] (that is computationally cheap to evaluate) to a function f will yield bounds that correspond to an expectation operator applied to an additive measure. Proposition 3. µ M is a sub-modular set function.

Proof. Following the ideas of [START_REF] Campos | Probability intervals: a tool for uncertain reasoning[END_REF] that prove super-modularity of µ M when ψ are probabilities, we will demonstrate the super-modularity of µ M , which by Proposition 1 is the same as proving sub-modularity of µ M . To do so, we consider arbitrary subsets A, B ⊆ X and we define the sets The idea of the proof is to define an additive set function κ ψ : 2 X → R and its mass function

C = A ∩ B and D = (A ∪ B) \ (A ∩ B). We have C ∩ D = ∅ and C ∪D = A∪B, as well as µ M (A∩B)+µ M (A∪B) = µ M (C)+µ M (C ∪D).
ψ such that µ M (C) = κ ψ (C) and µ M (C ∪D) = κ ψ (C ∪D) and ψ ≤ ψ ≤ ψ.
This is sufficient to prove sub-modularity of µ M , since

κ ψ (C) + κ ψ (C ∪ D) = κ ψ (B) + κ ψ (A) due to its additivity, and κ ψ (A) ≥ µ M (A) and κ ψ (B) ≥ µ M (B), follows from ψ ≤ ψ ≤ ψ.
There are four possible cases to study (two possible values for µ M (C) and two for µ M (C ∪ D)). All the cases are similar since the sets C and D are such that C ∩ D = ∅. We will therefore provide the details for only one of them.

Let us consider the specific case µ M (C) = x∈C ψ(x) and µ M (C ∪D) = κ(X )-

x ∈C∪D ψ(x). On each x ∈ X we are going to define the weights ψ(x) as follows:

-∀x ∈ C, ψ(x) = ψ(x), -∀x ∈ C ∪ D, ψ(x) = ψ(x),
This ensures that κ ψ coincide with µ M (C ∪ D) and µ M (C)

We are now going to see how we can assign masses d x := ψ(x) on elements of the subset D so that that ψ(x) ≤ ψ(x) = d x ≤ ψ(x) and x∈X = κ(X ). Observe that we still have to assign a total (potentially negative) mass κ(X )x ∈C∪D ψ(x) -x∈C ψ(x) on elements of D. To do so, consider first that, given the considered case, we have

µ M (C ∪ D) -µ M (C) = κ(X ) - x ∈C∪D ψ(x) - x∈C ψ(x). As, µ M (C) = x∈C ψ(x) ≥ κ(X ) -x ∈C ψ(x), we have µ M (C ∪ D) -µ M (C) ≤ κ(X ) - x ∈C∪D ψ(x) -κ(X ) + x ∈C ψ(x) ≤ x∈D ψ(x). Similarly, as µ M (C ∪ D) = κ(X ) -x ∈C∪D ψ(x) ≥ x∈C∪D ψ(x), we have µ M (C ∪ D) -µ M (C) ≥ x∈C∪D ψ(x) - x∈C ψ(x) = x∈D ψ(x). So x∈D ψ(x) ≤ µ M (C ∪ D) -µ M (C) ≤ x∈D ψ(x)
. This means that we can divide the difference µ M (C ∪ D) -µ M (C) (equal to the remaining mass to assign) on D and choose a real

d x for all x ∈ D such that ψ(x) ≤ d x ≤ ψ(x) and x∈D d x = µ M (C ∪ D) -µ M (C).
What Propositions 2 and 3 show is that previous results obtained for the probabilistic case (summative functions ψ summing up to one) readily extend to more general and complex cases, opening up new areas of applications such as the signal treatment one described in Section 4 and already suggested in Example 1.

Approximating sub-modular functions

We will now show that our approximation also works when one starts from a submodular function κ to approximate. There are at least two ways in which our proposed results can be useful: as an upper bound 4 , it can be useful in Branch and Bound techniques applied to sub-modular optimization [START_REF] Uematsu | An efficient branch-and-cut algorithm for submodular function maximization[END_REF]; as a trace on the singletons of a potentially unknown κ, it allows to build an approximation of κ while collecting only a 2n of its values.

We consider a sub-modular κ : 2 X → R. Recall that the conjugate set-function κ c of κ is ∀A, κ c (A) = κ(X ) -κ(A c ). It is easy to check that for any pair A, B ⊆ X , κ c is such that

κ c (A ∪ B) + κ c (A ∩ B) ≥ κ c (A) + κ c (B),
it is super-modular. κ is sub-modular so we have κ c ≤ κ. Now let us assume that we know only κ and κ c on each element x ∈ X . This is equivalent to know κ on each x and {x} c . In Example 2, this could amount to ask group member the value of every single food element, as well as the value of packages if we remove each food element, assuming that they aggregate values in a sub-modular way. By Equation ( 4), it means that adding a food element to an existing package is all the less interesting as the package grows bigger, which seems a reasonable assumption. For all x ∈ X let us define5 ψ κ (x) = κ c ({x}) and ψ κ (x) = κ({x}). The knowledge of ψ κ (x) and ψ κ (x) entails the following property:

Proposition 4. x∈X ψ κ (x) ≤ κ(X ) ≤ x∈X ψ κ (x). Proof. ψ k (x) ≤ ψ k (x) = κ({x}) = κ(X ) -κ c ({x} c ) ≤ κ(X ) -x =x ψ k (x ), so we have x∈X ψ k (x) ≤ κ(X ). ψ k (x) ≥ ψ k (x) = κ c ({x}) = κ(X ) -κ({x} c ) ≥ κ(X ) -x =x ψ k (x ), so we have κ(X ) ≤ x∈X ψ k (x).
As in the previous section, we now define two set functions

µ κ (A) = max{ x∈A ψ k (x), κ(X ) - x ∈A ψ k (x)}, (12)
µ κ (A) = min{ x∈A ψ k (x), κ(X ) - x ∈A ψ k (x)}, (13)
with the subscript κ denoting that they come from an initial sub-modular set function.

We will now proceed to show that µ κ and µ κ approximate κ and that they are respectively sub-and super-modular. Indeed, even if κ has a non-empty core C κ , we are not especially interested in positioning µ κ and µ κ with respect to this core, but with respect to the initial κ. Let us first show that µ κ and µ κ can act as lower/upper bounds.

Proposition 5. µ κ ≤ κ c ≤ µ κ .
Proof. Let us first prove that µ κ ≤ κ c . There are two possible values for µ κ (A):

-If µ κ (A) = x∈A ψ k (x) then µ κ (A) ≤ κ c (A) since κ is a sub-modular set func- tion. -If µ κ (A) = κ(X ) -x ∈A ψ k (x) then x ∈A ψ k (x) ≥ κ(A c ) = κ(X ) -κ c (A) which entails µ κ (A) ≤ κ(X ) -κ(X ) + κ c (A) = κ c (A).
The proof of κ c ≤ µ κ is similar.

Note that Proposition 1 also applies here, from which follows that µ κ and µ κ can also be used to approximate κ, as we can reproduce Proposition 5 proof for κ.

Corollary 1. µ κ ≤ κ c ≤ κ ≤ µ κ
We can now prove that µ κ enjoys the same interesting property as the function κ it approximates. Proposition 6. µ κ is a sub-modular set function i.e.,

∀A, B ⊆ X , µ κ (A ∪ B) + µ κ (A ∪ B) ≤ µ κ (A) + µ κ (B).
Proof. The proof is quite simple once we notice that C κ is not empty and that

ψ κ (x) = inf ψ∈Cκ ψ(x) ψ κ (x) = sup ψ∈Cκ ψ(x),
which can be obtained by applying Equations ( 9) and ( 7) to the indicator function I x (= 1 on x and zero everywhere else). Since the bounds ψ κ and ψ κ are atomic bounds over the set of additive measures C κ , one can apply Proposition 3.

Corollary 2. µ κ is super-modular.

Our results demonstrate that µ κ and µ κ are interesting approximations of a submodular function κ, possibly partially unknown. They are tight outer approximations, in the sense that as ψ κ (x) = κ c ({x}) and ψ κ (x) = κ({x}), any smaller intervals would not be outer approximations on the singletons. Finally, we can show an interesting link between the precision of the intervals [ψ κ , ψ κ ] and the functions µ κ , µ κ .

Proposition 7. For all A ⊆ X ,

|µ κ (A) -µ κ (A)| ≤ card(A)||u -l|| A,∞
where ||u -l|| A,∞ = max x∈A |ψ(x) -ψ(x)| and card(A)is the cardinality of A.

Proof. According to the definition of µ κ and µ κ , we have for possible cases.

µ κ (A) -µ κ (A) = -µ κ (A) = x∈A ψ k (x) ≤ κ(X ) -x ∈A ψ k (x) µ κ (A) = x∈A ψ k (x) ≥ κ(X ) -x ∈A ψ k (x). µ κ (A) -µ κ (A) = x∈A (ψ k (x) -ψ k (x)). -µ κ (A) = x∈A ψ k (x) ≤ κ(X ) -x ∈A ψ k (x) µ κ (A) = κ(X ) -x ∈A ψ k (x) ≥ x∈A ψ k (x). µ κ (A)-µ κ (A) = x∈A ψ k (x)-κ(X )+ x ∈A ψ k (x) ≤ x∈A (ψ k (x)-ψ k (x)). -µ κ (A) = κ(X ) -x ∈A ψ k (x) ≤ x∈A ψ k (x) µ κ (A) = x∈A ψ k (x) ≥ κ(X ) -x ∈A ψ k (x). µ κ (A)-µ κ (A) = κ(X )-x ∈A ψ k (x)-x∈A ψ k (x) ≤ x∈A (ψ k (x)-ψ k (x)). -µ κ (A) = κ(X ) -x ∈A ψ k (x) ≤ x∈A ψ k (x) µ κ (A) = κ(X ) -x ∈A ψ k (x) ≥ x∈A ψ k (x). µ κ (A)-µ κ (A) = κ(X )-x ∈A ψ k (x)-κ(X )+ x ∈A ψ k (x) ≤ x∈A (ψ k (x)- ψ k (x)).

So we have |µ

κ (A) -µ κ (A)| ≤ card(A)||u -l|| A,∞ .
Note that if κ is linear the previous result proves that µ κ = µ κ which is trivially true. More over, the further κ is from additivity ,the greater the difference |µ κ (A) -µ κ (A)| can be. This result also seems to connect our approximation with the total variation distance, routinely used in statistics [START_REF] Devroye | No empirical probability measure can converge in the total variation sense for all distributions[END_REF].

Finally, let us remark that in Example 2, we would have derived the same bounds if we assumed that the group computed the utility of a subset A according to a submodular aggregation function. In particular, this means that if we concentrate on atomic bounds only, it does not really matter whether µ and µ approximate a set M of additive set functions or a sub-modular set function (which itself can be considered as a nonempty set of dominated additive set functions), as both approximations will coincide if the bounds ψ, ψ are the same.

An application to signal convolution

We illustrate our approximation technique on an application in signal processing, as signal processing is strongly based on the notion of linear system, i.e. additive set functions.

Signal processing and additive set functions

We speak of continuous systems when the signals are continuous and of discrete systems when the signals are discrete. However, a signal is rarely discrete by nature. Discrete signals are generally continuous signals that have been discretized. For example, natural digital images are discretizations, by the retina of a camera, of continuous luminance images in real space. When they are naturally discrete, signals are often produced to mimic a discretized signal (synthetic images and sounds) with the aim of being used in the real world (special effects, electronic music, etc.).

Except in very rare cases, signal processing only really makes sense in continuous space. When the processing of the sound signal is used to recognize a sound for example (voice dictation) it is indeed recognition of the continuous signal that is involved via the recognition of the discrete signal from which it originates. So, doing digital signal processing generally means digitally performing a modification or analysis of the continuous signal from which it is being sampled. We limit ourselves here to linear signal processing.

The classical view consists in supposing that the processing of the signal, for which it is a question of finding a digital equivalent, is obtained thanks to a known continuous linear system (e.g. noise reduction filter (sound), texture modification (images), . . . ).

When a system (continuous or discrete) is linear, then the relation between input and output is completely defined by a continuous (ϕ) or discrete (ψ) signal called impulse response of the system and a linear operation called convolution. If the system is continuous, s the output of the system is written as a function of e the input:

s(t) = (e * ϕ)(t) = ∞ -∞ e(u).ϕ(t -u)du, (14)
and in the discrete case :

s(y) = (e * ψ)(y) = ∞ x=-∞ e(x).ψ(y -x). (15)
Most signal processing being inherently continuous, one needs to go from the continuous ideal impulse response ϕ(t), t ∈ R to a corresponding discrete impulse response ψ(x) x ∈ Z. One standard way to do that is to consider ψ(x) = β.ϕ(k.∆), ∆ being the sampling period, and β a normalization coefficient allowing to respect the dynamics of the system, i.e.

R ϕ(t)dt = x∈Z ψ(x). This approach is based on a modeling of the sampling by a Dirac impulse, which is very far from reality.

A more reasonable way to see the problem is to assume that the sampling was performed by a real system allowing to go from continuous space to discrete space whose impulse response η is of finite support and bound to the sampling period ∆. This is written ψ(x) = (ϕ * η)(x.∆) (see Figure [START_REF] Cohen | Cardinal extensions of the eu model based on the choquet integral[END_REF]). η is called the summative sampling kernel, summative meaning that it integrates to one ( R η(t)dt = 1) and its support is included in [-∆, ∆] (i.e. ∀t / ∈ [-∆, ∆], η(t) = 0). Note that the input signal e being generally finite, i.e. it is made of n samples, the expression (15) is generally written:

s(y) = (e * ψ)(y) = n x=1 e(x).ψ(y -x) = E κ ψy (e), (16)
with ∀y ∈ Z, ψ y (y) = ψ(y -x). When the operation of the equation ( 16) does not represent a system but a mathematical operation aimed at modifying a signal, ψ is more generally named convolution kernel.

Sampling and fuzzy transformation

A major problem posed by this approach is that the values of ψ strongly depend on the arbitrary choice of the chosen sampling kernel. For the discrete-continuous approach to be valid, the continuous convolution kernel must be digitized with the same sampling kernel as the discrete signal to be filtered. In most applications, this kernel is at best unknown or even, in the case of image processing, non-existent because the sampling is not invariant by translation: each pixel has its own impulse response.

To overcome this problem, Loquin et al. [START_REF] Loquin | On the granularity of summative kernels[END_REF] propose to replace the summative sampling kernel by a family of continuous sampling kernels (which are formally equivalent to probability distributions) represented by a possibility distribution π [START_REF] Dubois | When upper probabilities are possibility measures[END_REF], and called maxitive kernels. Within the framework of this theory, the triangular maxitive kernel plays a particular role since it allows to represent all the symmetric monomodal sampling kernels of lesser granularity [START_REF] Baudrit | Practical representations of incomplete probabilistic knowledge[END_REF].

This technique has been used in [START_REF] Strauss | Non-additive interval-valued f-transform[END_REF] under the name of imprecise fuzzy transformation to obtain an interval-valued impulse response [ψ(x), ψ(x)] (x ∈ Z) of all the values that could have been obtained by considering a family of sampling kernels represented by π.

We propose to extend the convolution operation defined in Equation ( 16) to the interval-valued impulse response [ψ, ψ] by:

[s(y), s(y)] = (e * [ψ, ψ])(y) = [E µ My (e), E µ My (e)], (17)
with µ My being defined by Equations [START_REF] Peter | Binary interactions and subset choice[END_REF], considering ∀x ∈ Z, ψ y (x) = ψ(y -x) and likewise for µ My .

Derivation of a discrete signal

The derivation of a discrete signal does not make sense. Indeed, only continuous signals can be differentiable. However, many signal analysis techniques are based on their derivation (extraction of contours or particular points in images, extraction of patterns in physiological signals, . . . ). In this context, it is not the discrete signal that we want to derive. What we want to do is to find an operation equivalent to a derivation of the continuous signal (before digitization) then sampling of the derivative obtained. It has been shown in [START_REF] Unser | B-spline signal processing. i. theory[END_REF] that this succession of operations is equivalent to convoluting the discrete signal with the sampling of the derivative of the kernel used for the reconstruction. One of the kernels often used to perform this operation is the Shen-Castan kernel [START_REF] Shen | Towards the unification of band-limited derivative operators for edge detection[END_REF] ∀t It is interesting to note that, in this application, µ My (X ) = µ My (X ) ≈ 0. Figure [START_REF] Baudrit | Practical representations of incomplete probabilistic knowledge[END_REF] shows the corresponding bounds, together with some discrete ψ obtained by the use of sampling kernels that are within the family represented by π, that is symmetric kernels with support [-δ, δ] (δ ∈]0, 1]). Such kernels includes, e.g., linear and quadratic kernels, as well as the Dirac kernel (whose result is in black in Figure [START_REF] Baudrit | Practical representations of incomplete probabilistic knowledge[END_REF]).

∈ R, φ(t) = - ln(β) 2 .β |t| , (18)

Conclusion

In this paper, we have discussed how having only atomic bounds can be used to efficiently approximate very generic sets of very generic additive or sub-modular measures (i.e., non-monotonic with arbitrary bounds).

We have shown how they can successfully be applied to perform robust signal processing. Our next step will be to investigate how our results regarding sub-modular and sec. set functions can be instrumental to represent non-monotonic rankings of sets or as approximation tools in optimisation problems using sub-modular functions.

  So we are going to prove the property considering C and D subsets such that C ∩ D = ∅.

Fig. 1 :

 1 Fig. 1: Sampling of the impulse response ϕ with the kernel η.

where β ∈]0, 1 [ 2 ,

 12 is a form factor of the kernel. The derivation of the kernel φ leads to the continuous kernel ∀t ∈ R, ϕ(t) = -sign(t). ln(β)2 2 .β |t|Let us assume, without any loss of generality, that the sampling step ∆ = 1. Then, sampling ϕ as proposed in Section 4.2 by using a possibility distribution π leads to the interval-valued impulse response of the form [ψ, ψ]:if x > 0, ψ(x) = -ln(β).(β |x+1| -β |x| ) 2 , ψ(x) = -ln(β).(β |x| -β |x-1| ) and if x < 0, ψ(x) = ln(β).(β |x| -β |x-1| ) 2 , ψ(x) = ln(β).(β |x| -β |x| ) 2 .

3 Fig. 2 :

 32 Fig. 2: ψ (in blue), ψ (in red), 8 impulse responses for different kernels within the family represented by π (in magenta), and the response of the Dirac kernel (in black).

Figure ( 3 )

 3 Figure[START_REF] Benabbou | Interactive optimization of submodular functions under matroid constraints[END_REF] and Figure[START_REF] Prade | Possibility Theory. An Approach to Computerized Processing of Uncertainty[END_REF] respectively show a time-attenuated signal and its derivative bounds obtained by our approach, together with some sampled filters. We can see on Figure (4) that our approach makes a good job, providing wide bounds when variations are important (bottom box in Figure (4)), while providing very tight bounds when variation decreases along time (upper box).

  sec.

Fig. 3 :

 3 Fig. 3: Synthesis signal composed of a sum of attenuated sinusoids.

Fig. 4 :

 4 Fig. 4: Original signal derived by the imprecise Shen filter (lower bound in red, upper bound in blue) superimposed with derived signals obtained by various sampled Shen filters.

Table 1 :

 1 Table of scores of Example 2

or lower bound if we consider the complementary measure

The reason for also calling them ψ is clarified in the proof of Proposition