**1. Introduction**

The oceans abound with natural physical sounds (from wind, rain, polar ice, and seismic activity), biological sounds (from crustaceans, fishes, and marine mammals), and anthropogenic sounds (from transport, construction, offshore exploration, and mining). Soundscapes naturally change over time because of temporal cycles in weather (e.g., cyclones and annual monsoon [1,2]) and animal behaviour (e.g., diurnal foraging patterns, lunar spawning, seasonal mating, and annual migration [3–6]). However, in many habitats, soundscapes further change with patterns of human presence (e.g., temporary construction or summer recreation [7]) and some have changed steadily over time with increasing intensity of anthropophony (e.g., due to shipping [8]).

In 1996, the European Commission identified air-borne noise as one of the main terrestrial environmental issues in Europe, having been neglected compared to chemical pollution [9]. Subsequently, the Commission enacted sound mapping as an important step to assess and manage sound exposure levels in urban areas [10]. A little later, the issue of underwater ocean noise received similar attention, being declared a pollutant, and with underwater sound monitoring and mapping being suggested [11]. Nowadays, underwater noise footprints of individual anthropogenic operations are commonly mapped for environmental impact assessments (e.g., [12–14]). Longer-term, large-scale marine

**Citation:** Erbe, C.; Schoeman, R.P.; Peel, D.; Smith, J.N. It Often Howls More than It Chugs: Wind versus Ship Noise Under Water in Australia's Maritime Regions. *J. Mar. Sci. Eng.* **2021**, *9*, 472. https://doi.org/ 10.3390/jmse9050472

Academic Editor: Michele Viviani

Received: 30 March 2021 Accepted: 25 April 2021 Published: 27 April 2021

**Publisher's Note:** MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

**Copyright:** © 2021 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (https:// creativecommons.org/licenses/by/ 4.0/).

sound mapping has focussed on ship noise [15–17], but may also include other sources, such as seismic airguns and explosives [18].

Shipping is a global contributor to ocean noise and, over the past five decades, has caused a steady increase in underwater low-frequency (10–100 Hz) ambient sound levels in many marine regions [19–24]. This is of concern, because ship noise causes behavioural and acoustic responses, auditory masking, and stress in marine animals [25–28]. Hence, various studies have mapped ship noise and overlain the resulting maps with marine habitat maps to identify areas of concern (hotspots; high animal density and high noise) [29] and areas of opportunity (high animal density and low noise) [30] for marine spatial planning.

A problem with ship noise maps is that they often lack validation against in situ measurements. These maps may have several sources of error in the ship positions and routes, source spectra and levels, sound propagation models, and hydro- and geoacoustic parameters required by the models. As well, the spatial (depth and range) and temporal grid over which the models operate introduces uncertainty. In fact, lack of knowledge on the physical environment (i.e., hydroacoustic parameters of the water and geoacoustic parameters of the upper seafloor) is often the limiting factor in sound propagation model accuracy [14,31]. Model validation is essential to confirm accuracy and to support the use of a sound map for management decisions [32].

Finally, underwater anthropogenic noise needs to be put into context. How does it compare to natural, pervasive noise as from wind? Sertlek et al. [18] found that shipping inserted the greatest amount of acoustic energy into the Dutch North Sea and far exceeded that of wind. Similarly, Farcas et al. [32] showed that ship noise exceeded wind noise under water near major ports and shipping lanes, and around industrial sites in the Northeast Atlantic. However, southern hemisphere oceans have a reputation of being less impacted by anthropogenic sounds, largely due to a lower ship density [33]. Thus, wind may supersede ship noise in parts of the southern oceans. Here, we model underwater sound in the Australian Exclusive Economic Zone (EEZ) from ships and wind over a 6-month austral winter period, and validate the model with 6-month recordings from 25 stations. We chose the winter months as this is the peak of baleen whale presence (e.g., [34–36]). The aim is to enable a better understanding of where shipping noise is likely to be a significant contributor to the marine soundscape and thus, a potential stressor to marine life.

#### **2. Materials and Methods**

In a nutshell, we used ship tracks from Automatic Identification System (AIS) logs and underwater source spectra from the literature. On a 5 km × 5 km grid over the Australian EEZ, we modelled underwater sound propagation from all source cells (i.e., cells that contained ships) to all surrounding receiver cells over a 100 km radius. We then integrated sound exposure over the austral winter. The computational effort was managed by (1) splitting the EEZ into 28 acoustic zones, in which sound propagates in similar ways [37], hence, where a similar model may be set up, and (2) using a neural network to cluster all source-receiver transects within a zone into 64 groups of bathymetry transects, and modelling sound propagation only for cluster centroids. An overview of the process step-by-step is given in Appendix A. Wind noise was not propagated, but simply computed based on hourly wind speed data in each cell.

All GIS analysis was done with a combination of ArcMap (version 10.5, ESRI, Redlands, CA, USA) and R (Version 4.03, R Core Team, Vienna, Austria). Noise modelling and validation was done in MATLAB (Version 2020b, The MathWorks Inc., Natick, MA, USA). We commenced with a GIS layer of the Australian marine bathymetry, gridded to 5 km × 5 km [38].

#### *2.1. Ship Data*

Data on ship type, size, position, and speed were obtained from Automatic Identification System (AIS) logs managed by the Australian Maritime Safety Authority (AMSA). AIS data were extracted for the period 1 April 2015–30 September 2015. Ships were grouped

into five classes based on their length only (and not by type or function; e.g., tanker versus passenger ferry): <25 m, ≥25–<50 m, ≥50–<100 m, ≥100–<200 m, and ≥200 m. The regularity of the AIS location reporting depends on location and time, and the data we used was subsampled to provide locations, at most, every 5 min. From this data, ship tracks were interpolated by dead reckoning to intervals of 60 s if two successive AIS positions met criteria based on the time between the polls and the straightness of the vessel's path (as per Appendix B in [39]). For each ship track, the time spent in each grid cell was computed, and time was summed over all ships belonging to the same class, yielding a grid of cumulative time spent in each cell over the 6-month period, by ship class.

Underwater ship noise source spectra were taken from the Research Ambient Noise Directionality (RANDI) model [40] and integrated into full-octave bands: ≥10–<20 Hz, ≥20–<40 Hz, ≥40–<80 Hz, ≥80–<160 Hz, ≥160–<320 Hz, ≥320–<640 Hz, ≥640–<1280 Hz, and ≥1280–<2560 Hz (Figure 1a). The broadband (10 Hz–2.6 kHz) source levels were: 148, 160, 172, 187, and 193 dB re 1 μPa m, for the five classes, respectively.

**Figure 1.** (**a**) Underwater ship noise source levels as full-octave band levels (OBL) for ships of lengths <25 m (Class 1), ≥25–<50 m (Class 2), ≥50–<100 m (Class 3), ≥100–<200 m (Class 4), and ≥200 m (Class 5); (**b**) Power spectral density levels (PSD) of wind noise under water at wind speeds 1–3 kn (Beaufort 1; curve 1), 4–6 kn (Beaufort 2; curve 2), 7–10 kn (Beaufort 3; curve 3), 11–21 kn (Beaufort 4–5; curve 4), 22–47 kn (Beaufort 6–9; curve 5), and ≥48 kn (≥Beaufort 10; curve 6).

#### *2.2. Wind Data*

Hourly data on surface wind speed (10 m altitude) were obtained over a similar 6-month period (1 April 2012–30 September 2012) from the Bureau of Meteorology and CSIRO [41], based on the NCEP Climate Forecast System [42]. The data varied in spatial resolution (4, 10, and 24 arcminute grids), which we projected and re-sampled to a 5 km × 5 km UTM grid. Over all grid cells, wind speed varied between 0.5 and 30 m/s (i.e., 1–58 kn). Wind speeds were binned to match the sea states represented in the 'Wenz curves' [43] and noise spectra were assigned to each wind speed bin (Figure 1b). The 'Wenz curves' were converted to linear power spectral density, then integrated over frequency, before applying 10log10 to yield broadband mean-square sound pressure levels. Expressed as root-mean-square sound pressure levels, the associations became: 79 dB re 1 μPa for ≥1–<4 kn, 87 dB re 1 μPa for ≥4–<7 kn, 92 dB re 1 μPa for ≥7–<11 kn, 99 dB re 1 μPa for ≥11–<22 kn, 105 dB re 1 μPa for ≥22–<48 kn, and 113 dB re 1 μPa for ≥48 kn.

#### *2.3. Acoustic Zones*

The Australian EEZ had previously been broken up into 28 'acoustic zones' (Figure 2), whereby each zone was characterised by a unique set of hydroacoustic parameters of the water (i.e., sound speed profiles), geoacoustic parameters of the seafloor (i.e., thickness of

the sediment layer, density, compressional sound speed and attenuation, and shear sound speed and attenuation), and bathymetric parameters (i.e., water depth and slope) [37]. The idea was to set up one sound propagation environment for each zone and then model all of the ships in that zone.

**Figure 2.** Map of the marine acoustic zones of the Australian EEZ [37].

#### *2.4. Source-Receiver Transects*

Within each zone, all of the grid cells that had a cumulative ship time > 0 (no matter the ship class) were identified. From each of these 'source cells', 36 radials were cast at 10-degree intervals. The bathymetry was extracted along each of these radials in 5 km steps out to a maximum range of 100 km, yielding thirty-six 100 km source-receiver transects around each source. A bathymetry reading at 2.6 km range from the centre of each source cell was added, representing the mean distance between two random points inside a square (i.e., 0.5214 times the edge length [44]). If a transect hit land, all subsequent bathymetry samples were set to 'not a number' along this 100 km transect. If a source sat near a zone boundary, then the 100 km transects were extracted with bathymetry from the neighbouring zone or from a 100 km buffer around the outside of the EEZ.

All of the transects from all of the source cells (of all ship classes) within a zone were passed to an unsupervised Kohonen neural network (i.e., a self-organizing map, SOM) with 900 neurons [45] (also see [13], where this SOM was previously used to cluster bathymetry transects). The neural net sorted the transects into 900 groups, based on their bathymetric shape. Further grouping was achieved by k-means clustering allowing for 64 clusters [46]. Cluster centroids were computed as the arithmetic mean of all transects within one cluster (see examples for Zone 16 in Figure 3). Sound propagation was modelled along each of the 64 centroids within one zone.

**Figure 3.** Graphs of all 98,532 source-receiver transects of Zone 16 plotted by cluster (grey), with centroid bathymetries shown in black. X-axes are range (km) and Y-axes are depth below the sea surface (m). Note the changing Y-ranges.

#### *2.5. Sound Propagation Model*

Sound propagation over each centroid bathymetry was modelled with RAMGeo in AcTUP V2.8 [47] (https://cmst.curtin.edu.au/products/underwater/ accessed on: 27 March 2021) based on zone-specific acoustic environments consisting of three layers: the water column, an unconsolidated surface sediment layer, and a consolidated calcarenite sediment layer as a half space. Water column parameters included the zone's mean sound speed profile [37] and water density profile. Representative temperature and salinity data were extracted from the World Ocean Atlas [48,49] to calculate water densities based on the UNESCO formula for sea water density [50]. Unconsolidated surface sediment layers throughout the EEZ comprised predominantly fine material (silt-sand) with sufficiently low shear wave speeds (<250 m/s) to allow modelling as a fluid. Hence, unconsolidated surface sediment parameters only included the zone-specific layer thickness, compressional sound speed, compressional wave attenuation, and density [37]. Surface sediment layer thickness was estimated as 0.5 m for zones within the sediment-starved carbonate platform [51]. Surface sediment thickness in the remaining zones appears variable [52–57], and so was modelled as 2 m.

In contrast to the surface sediment layer, calcarenite acts as an elastic material with a shear sound speed of 1400 m/s resulting in important propagation effects [51]. While RAMGeo is a parabolic equation for fluid seabeds, reasonable results can be obtained by using an equivalent fluid approximation with reflection coefficients representative of the elastic model [58]. The procedure to find an equivalent fluid approximation for each zone included (a) creating an environment with calcarenite as an elastic material (see [51] for geoacoustic properties), (b) creating an environment with calcarenite as a fluid layer starting with a compressional sound speed of 1250 m/s and a compressional wave attenuation of 4.5 dB/λc, (c) calculating the reflection coefficients for each modelled frequency for both environments with the reflection coefficient model BOUNCE [59], and (d) adjusting the fluid layer parameters until a representative equivalent had been reached (Figure 4).

**Figure 4.** Example of the reflection coefficients of an elastic model (blue line) and a representative fluid model (green line) calculated with BOUNCE at 125 Hz.

RAMGeo modelled sound propagation loss for the centre frequencies of eight fulloctave bands between 10 and 2000 Hz over a range of 100 km and up to 7100 m depth, which is more than the maximum water depth (6388 m) of the EEZ. The depth and range resolutions were 10 m. The source depth for all ships was chosen as 5 m below the sea surface. An example RAMGeo output at 640 Hz for the 64 centroid bathymetries of Zone 16 is shown in Figure 5. The bathymetry itself is just visible as a black line, below which propagation loss was greatest. The colours vary from 60 dB propagation loss (dark red) to 110 dB (dark blue). Several patterns are obvious: Convergence zones appeared over all deep bathymetries, leading to low propagation loss (i.e., high received levels) near the sea surface about 60 km from the source (i.e., clusters 6, 11, 26, 32, 41, 52, 53, 55, 59, 60, and 64). Over upwards-sloping bathymetries, propagation loss was greatest (e.g., clusters 1, 4, 19, 35, 37, 54, and 58). Over downwards-sloping bathymetries, sound may reflect into the deep-ocean sound channel, which has an axis at about 1 km depth off Australia. Once inside the channel, sound may propagate over vast ranges at very little additional loss because of no further seafloor (and to a lesser extent, sea surface) interactions (i.e., clusters 12, 13, 25, 36, and 49). Finally, RAMGeo does not include frequency-dependent absorption and so this additional loss was applied outside of and after RAMGeo [60].

**Figure 5.** Plots of propagation loss (PL) as a function of range and depth along the 64 centroid bathymetries from Zone 16, for a frequency of 640 Hz. The darkest red corresponds to 60 dB and the darkest blue to 110 dB. X-axes are range (km) and Y-axes are depth below the sea surface (km); both are scaled linearly.

### *2.6. Accumulation of Received Levels*

Within each zone, one ship class was treated at a time. The source cells corresponding to one ship class were found, thirty-six 100 km radials were cast at 10-degree intervals around each source cell, and bathymetry was extracted along each radial and sampled in 5 km steps; the mean distance between a ship and a receiver of 2.6 km within the source cell was inserted at the beginning. In other words, source cells were assigned a received level at 2.6 km. Then, stepping through the source cells for this ship class in this zone, for each of the 36 radial transects, the best matching centroid bathymetry was found. In fact, as the SOM had been trained with all source-receiver transects from all ship classes in this zone, finding the 'best matching centroid' reduced to simply looking up into which cluster this transect originally went. Then, for each frequency modelled, propagation loss (*PL*) along this centroid was recovered and subtracted from the corresponding octave band

source level (*SL*) plus the cumulative time in dB re 1 s (with *T* in units of second) of this ship class in this source cell, yielding received levels (*RL*):

$$RL = SL + 10\log\_{10}(T) - PL$$

In this equation, *SL* is a number, the octave band source level expressed as a meansquare pressure level [dB re 1 μPa2] at the modelled frequency. The duration term is also a number [dB re 1 s]. *PL* [dB], however, is a matrix with values in depth and range. Therefore, *RL* is a matrix containing received sound exposure levels (*SEL*) [dB re 1 μPa2s] as a function of range and depth. There was one such *RL* matrix for each frequency.

To change from the polar coordinates (in which sound propagation was modelled) to the Cartesian grid of the EEZ, *RL* was interpolated to the 5 km grid of the EEZ, at each depth. Including frequency as an additional dimension, this yielded a 4-dimensional matrix of longitude, latitude, depth, and frequency, covering the entire EEZ. The matrix was populated cumulatively by summing sound exposure (i.e., in linear, not logarithmic terms) over all 36 transects about each source cell, and then over all source cells, before taking 10 log10 again to yield cumulative sound exposure levels (*C-SEL*).

#### *2.7. Ship Noise Map*

Broadband sound exposure levels were computed by summing sound exposure over frequency, thereby reducing the matrix to three dimensions, then converting to dB. A further reduction to two dimensions was achieved by finding the maximum sound exposure level over the top 200 m, representing the 'worst case' of exposure for animals that dive over this depth [61]. One such map is presented for each ship class, as well as cumulatively over all five classes. These maps of cumulative sound exposure level were accumulated over six months encompassing the austral winter. They can be read as average mean-square sound pressure level (*SPL*) maps by subtracting the 6-month duration in dB re 1 s:
