\documentclass[10pt]{article}

\usepackage{fullpage}
\usepackage{setspace}
\usepackage{parskip}
\usepackage{titlesec}
\usepackage[section]{placeins}
\usepackage{xcolor}
\usepackage{breakcites}
\usepackage{lineno}
\usepackage{hyphenat}





\PassOptionsToPackage{hyphens}{url}
\usepackage[colorlinks = true,
            linkcolor = blue,
            urlcolor  = blue,
            citecolor = blue,
            anchorcolor = blue]{hyperref}
\usepackage{etoolbox}
\makeatletter
\patchcmd\@combinedblfloats{\box\@outputbox}{\unvbox\@outputbox}{}{%
  \errmessage{\noexpand\@combinedblfloats could not be patched}%
}%
\makeatother


\usepackage{natbib}




\renewenvironment{abstract}
  {{\bfseries\noindent{\abstractname}\par\nobreak}\footnotesize}
  {\bigskip}

\titlespacing{\section}{0pt}{*3}{*1}
\titlespacing{\subsection}{0pt}{*2}{*0.5}
\titlespacing{\subsubsection}{0pt}{*1.5}{0pt}


\usepackage{authblk}


\usepackage{graphicx}
\usepackage[space]{grffile}
\usepackage{latexsym}
\usepackage{textcomp}
\usepackage{longtable}
\usepackage{tabulary}
\usepackage{booktabs,array,multirow}
\usepackage{amsfonts,amsmath,amssymb}
\providecommand\citet{\cite}
\providecommand\citep{\cite}
\providecommand\citealt{\cite}
% You can conditionalize code for latexml or normal latex using this.
\newif\iflatexml\latexmlfalse
\providecommand{\tightlist}{\setlength{\itemsep}{0pt}\setlength{\parskip}{0pt}}%

\AtBeginDocument{\DeclareGraphicsExtensions{.pdf,.PDF,.eps,.EPS,.png,.PNG,.tif,.TIF,.jpg,.JPG,.jpeg,.JPEG}}

\usepackage[utf8]{inputenc}
\usepackage[ngerman,greek,english]{babel}



\usepackage{float}








\begin{document}

\title{Combined genomics to discover genes associated with tolerance to soil
carbonate}



\author[1]{Charlotte Poschenrieder}%
\author[1]{Silvia Busoms}%
\author[1]{Laura Pérez-Martín}%
\author[1]{Joana Gelabert}%
\author[2]{Xin-Yuan Huang}%
\author[3]{Levi Yant}%
\author[1]{Roser Tolrà}%
\author[4]{David Salt E}%
\affil[1]{Universitat Autonoma de Barcelona Departament de Biologia Animal de Biologia Vegetal i d'Ecologia}%
\affil[2]{Nanjing Agricultural University State Key Laboratory of Crop Genetics and Germplasm Enhancement}%
\affil[3]{University of Nottingham School of Life Sciences}%
\affil[4]{University of Nottingham School of Biosciences}%


\vspace{-1em}



  
  \date{May 18, 2024}


\begingroup
\let\center\flushleft
\let\endcenter\endflushleft
\maketitle
\endgroup





\selectlanguage{english}
\begin{abstract}
Carbonate-rich soils limit plant performance and crop production.
Previously, local adaptation to carbonated soils was detected in wild
\emph{Arabidopsis thaliana} accessions, allowing the selection of two
demes with contrasting phenotypes: A1 (carbonate tolerant, c+) and T6
(carbonate sensitive, c-). Here, A1 \textsubscript{(c+)} and T6
\textsubscript{(c-)} seedlings were grown hydroponically under control
(pH 5.9) and bicarbonate conditions (10 mM NaHCO \textsubscript{3}, pH
8.3) to obtain ionomic profiles and conduct transcriptomic analysis. In
parallel, A1 \textsubscript{(c+)} and T6 \textsubscript{(c-)} parental
lines and their progeny were cultivated on carbonated soil to evaluate
fitness and segregation patterns. To understand the genetic architecture
beyond the contrasted phenotypes a bulk segregant analysis sequencing
(BSA-Seq) was performed. Transcriptomics revealed 208 root and 2503 leaf
differentially expressed genes (DEGs) in A1 \textsubscript{(c+)}
\emph{vs} T6 \textsubscript{(c-)} comparison under bicarbonate stress,
mainly involved in iron, nitrogen and carbon metabolism, hormones, and
glycosylates biosynthesis. Based on A1 \textsubscript{(c+)} and T6
\textsubscript{(c-)} genome contrasts and BSA-Seq analysis, 69 genes
were associated with carbonate tolerance. Comparative analysis of
genomics and transcriptomics discovered a final set of 18 genes involved
in bicarbonate stress responses that may have relevant roles in soil
carbonate tolerance.%
\end{abstract}%



\sloppy


Combined genomics to discover genes associated with tolerance to soil
carbonate

Silvia Busoms\textsuperscript{1}, Laura P\selectlanguage{ngerman}érez-Martín\textsuperscript{1},
Joana Terés\textsuperscript{1}, Xin-Yuan Huang\textsuperscript{2}, Levi
Yant\textsuperscript{3}, Roser Tolrà\textsuperscript{1}, David E
Salt\textsuperscript{4}, Charlotte Poschenrieder\textsuperscript{1*}.

\textsuperscript{1} Plant Physiology Laboratory, Bioscience Faculty,
Universitat Autònoma de Barcelona. C/ de la Vall Moronta s/n, E-08193
Bellaterra, Barcelona (Spain)\textsuperscript{2} State Key Laboratory of
Crop Genetics and Germplasm Enhancement, College of Resources and
Environmental Sciences, Nanjing Agricultural University, Nanjing 210095,
China\textsuperscript{3} Future Food Beacon of Excellence \& School of
Life Sciences, University of Nottingham, Nottingham, NG7 2RD,
UK.\textsuperscript{4}Future Food Beacon of Excellence \& School of
Biosciences, University of Nottingham, Sutton Bonington LE12 5RD, UK.

* Corresponding author: Charlotte Poschenrieder

\textbf{*} Correspondence:\emph{Charlotte.Poschenrieder@uab.es}

\textbf{Abstract}

Carbonate-rich soils limit plant performance and crop production.
Previously, local adaptation to carbonated soils was detected in
wild\emph{Arabidopsis thaliana} accessions, allowing the selection of
two demes with contrasting phenotypes: A1 (carbonate tolerant, c+) and
T6 (carbonate sensitive, c-). Here, A1\textsubscript{(c+)} and
T6\textsubscript{(c-)} seedlings were grown hydroponically under control
(pH 5.9) and bicarbonate conditions (10 mM NaHCO\textsubscript{3}, pH
8.3) to obtain ionomic profiles and conduct transcriptomic analysis. In
parallel, A1\textsubscript{(c+)} and T6\textsubscript{(c-)} parental
lines and their progeny were cultivated on carbonated soil to evaluate
fitness and segregation patterns. To understand the genetic architecture
beyond the contrasted phenotypes a bulk segregant analysis sequencing
(BSA-Seq) was performed. Transcriptomics revealed 208 root and 2503 leaf
differentially expressed genes (DEGs) in A1\textsubscript{(c+)}\emph{vs}
T6\textsubscript{(c-)} comparison under bicarbonate stress, mainly
involved in iron, nitrogen and carbon metabolism, hormones, and
glycosylates biosynthesis. Based on A1\textsubscript{(c+)} and
T6\textsubscript{(c-)} genome contrasts and BSA-Seq analysis, 69 genes
were associated with carbonate tolerance. Comparative analysis of
genomics and transcriptomics discovered a final set of 18 genes involved
in bicarbonate stress responses that may have relevant roles in soil
carbonate tolerance.

\textbf{Keywords:} Arabidopsis, soil carbonate, bicarbonate stress,
BSA-Seq; transcriptomics.

Introduction

Soil alkalinity is a highly stressful environmental factor that limits
plant growth and crop yield {[}1-2{]}. High pH soils are present in 30\%
of the Earth surface specially located in areas with arid and semiarid
climate. The pH values \selectlanguage{english}\selectlanguage{english}of most calcareous soils are within the range
of 7.5 to 8.5 {[}3{]} with bicarbonate concentrations between 5-35 mmol
L\textsuperscript{-1}. The main anions present in excess in alkaline
soils are HCO\textsubscript{3}\textsuperscript{-} and
CO\textsubscript{3}\textsuperscript{2-} {[}4{]}. Soil carbonates act as
pH buffers and play an important role in rhizosphere processes,
hampering nutrient availability to plants. The low availability of
nitrogen (N), phosphorus (P) and micronutrients such as iron (Fe), zinc
(Zn), and manganese (Mn) produce nutrient deficiencies in many plant
species cultivated on carbonated soils {[}5-7{]}. Moreover, the high pH
surrounding plant roots can alter the membrane potential {[}8{]} and
inhibit the uptake of essential ions, and thus further contribute to
nutrient deficiency in sensitive plants.

The uptake of bicarbonate can cause important metabolic disturbances due
to both alteration of cell pH homeostasis and dark fixation of inorganic
carbon. Bicarbonate can be quickly incorporated into organic acids,
mainly malate, causing inhibition of mineral nutrient transport from
roots to shoot, reduction of root growth, and oxidative stress in
sensitive calcifuges, but not in tolerant calcicole species {[}9,10{]}.
Due to these multiple direct and indirect injuries, sensitive species
exhibit a complex syndrome of stress symptoms, including
morpho-anatomical changes in roots, disturbance of water relations,
signs of iron and/or zinc deficiency, reduction of total photosynthetic
pigments and photosynthetic activity; accumulation of osmo-protectants,
soluble sugars, and organic acids; and activation of the biosynthesis of
antioxidant enzymes {[}11,12,13{]}.

A better understanding of plant alkaline tolerance mechanisms and
cultivation of new varieties of alkali-tolerant crops is needed to
improve carbonated soils and increase food production {[}14{]}. However,
current knowledge of the bicarbonate stress response of plants is
limited. Most studies have been conducted on roots of crop species such
as \emph{Glycine max} or \emph{Oryza sativa} which are species with
moderate bicarbonate tolerance {[}15, 16{]}. Nonetheless, after several
decades of effort still a proper model is missing to understand the
adaptive mechanism to carbonate stress and the tolerance mechanism at
the molecular-genetic level. Early experiments revealed that organic
acid accumulation in response to bicarbonate occurred both in sensitive
and tolerant species and it was speculated that differential
compartmentation of organic acids may play a role in the tolerance
{[}10{]}. In fact, several anion transporters, and gene families such as
ALMT, NRT/POT and SLAHs were related to bicarbonate response
in\emph{Glycine max} {[}16,17{]}.

During recent years, combined omics are applied to successfully
characterize plant tolerance responses to complex stress factors (e.g.
{[}18{]}). Resequencing-based Bulk Segregation Analysis (BSA) utilizes a
strategy of pooling individuals with extreme phenotypes to conduct
rapidly linked marker screening {[}19{]}. Whole-genome studies are
useful to detect signature of selection, QTLs, and SNP markers. However,
some important changes affecting the regulation of key genes can remain
hidden. The combination of genomic and transcriptomic studies allows a
more accurate screening of the candidate genes involved in specific
physiological processes. Here, we took advantage of the natural
variation present in two \emph{A. thaliana} populations with contrasting
phenotypes of soil alkalinity tolerance {[}5{]} to highlight the loci
involved in bicarbonate stress responses and adaptation to carbonated
soils by combining BSA-Seq and RNA-Seq technologies.

Materials and Methods

\emph{Plant material}

In previous studies, natural variation of \emph{Arabidopsis
thaliana}populations from Catalonia {[}20{]} were tested in a multi-year
small scale common garden under carbonated and non-carbonated soils
{[}5{]}. Seeds of 6 extreme phenotype plants from
A1\textsubscript{(c+)}, a soil carbonate tolerant deme, and 6 extreme
phenotype plants from T6\textsubscript{(c-)}, a soil carbonate sensitive
deme, were collected from the last reciprocal transplant experiment
performed in 2015 and stored under fresh (4o C) and dry conditions until
the beginning of hydroponic and soil cultivation experiments. A1xT6 and
T6xA1 crosses were generated using the 6 A1\textsubscript{(c+)} and 6
T6\textsubscript{(c-)} lines (generating 12 crossing lines) and
self-crossed for three generations (F1, F2, F3). Col-0 seeds and seeds
of 18 T-DNA mutants with Col-0 background were purchased from NASC/ABRC
{[}21{]}. Segregating lines were genotyped using specific primers and
multiplied for T2 generations. The non-expression of the genes was
validated through qPCR (Dataset\_SD1).

\emph{Growth conditions~}

Selected seeds were surface sterilized by soaking in 70\% (v/v) ethanol
for 1 min, suspended in 30\% (v/v) commercial Clorox bleach and 1 drop
of Tween-20 for 5 min and rinsed 5 times in sterile 18 M\selectlanguage{greek}Ω\selectlanguage{english}milli-Q
water.~Seeds were kept in water, at 4 \selectlanguage{ngerman}ºC, in dark conditions for 48h to
synchronize germination.

\emph{Hydroponic culture} : \emph{A. thaliana} seeds from 6
A1\textsubscript{(c+)} and 6 T6\textsubscript{(c-)} lines were sown in
0.2 ml tubes containing 0.6 \% phyto agar (Duchefa, Harlem, The
Netherlands) prepared in nutrient solution 1/2 Hoagland solution (HS),
pH 5.9 in a growth chamber (12 light/12 dark hours, 150 \selectlanguage{greek}µ\selectlanguage{english}mol
cm\textsuperscript{-2}\selectlanguage{ngerman}·s\selectlanguage{english}\textsuperscript{-1}, 40\% humidity and 25deg C).~
The bottom of the tubes containing seedlings was cut off and the tubes
were placed in 150 ml hydroponic container with aerated nutrient
solution 1/2 Hoagland, pH 5.9. After 15 days post-germination, treatment
was applied, and seedlings were separated in different sets. The
treatments consisted of control (1/2 HS at pH 5.9), high pH (1/2 HS at pH
8.3), and bicarbonate (1/2HS at pH 8.3 and 10 mM NaHCO\textsubscript{3}).
Solutions were buffered with different proportions of MES and BTP
(BIS-TRIS propane) depending on the final pH.~ Plants of
A1\textsubscript{(c+)}, T6\textsubscript{(c-)} used for transcriptomic
analysis were collected after 48 hours under treatment. Plants of
A1\textsubscript{(c+)}, T6\textsubscript{(c-)}, Col-0 and the 18 T-DNA
knockout mutants were collected after 10 days under treatment. Image
analysis for root length measurements were performed with Image J,
software16 {[}22{]}.

\emph{Greenhouse experiments} : four cultivations were conducted
consecutively under the same conditions: (1) 6 A1\textsubscript{(c+)}, 6
T6\textsubscript{(c-) ,} 6 A1xT6 and 6 T6xA1 F1 lines; (2) 5
A1\textsubscript{(c+)}, 2 T6\textsubscript{(c-) ,} 10 A1xT6 and 4 T6xA1
F2 families; (3) 5 A1\textsubscript{(c+)}, 2 T6\textsubscript{(c-),} 66
A1xT6 and 28 T6xA1 F3 families; (4) Col-0 and the 18 T-DNA mutants of
candidate genes. Plants were cultivated in carbonated soil (LP) mixed
with perlite (3:1) (Table 1) in 6x6cm square pots under semi-controlled
conditions (12-h light/12-h dark photoperiod, temperature between 15-25
oC). Five seeds of each line were sown in pots and distributed randomly
in the greenhouse. Two weeks after germination, seedlings were thinned
out so that only one plant per pot was left. Irrigation with distilled
water was applied twice a week at the bottom of the trays. For the third
cultivation, 2 leaves per plant were collected 20 days after germination
and stored at -80 \textsuperscript{o}C for subsequent DNA extraction and
sequencing of the selected plants. Rosette diameter and bolting time
were monitored, and silique number was counted at maturity in all the
experiments.

\emph{Ionomic analysis}

Root ionome was assessed in plants submitted for 10 days to treatment
conditions.~Roots were carefully washed with 18 M\selectlanguage{english} to remove nutrients
from solution. Plant material was dried for 4 days at 60 oC.
Approximately 0.1g was used to perform an open-air digestion in Pyrex
tubes using 0.7 mL concentrated HNO\textsubscript{3}~at 110 oC for 5 h
in a hot-block digestion system (SC154-54-Well Hot Block\selectlanguage{english}, Environmental
Express, SC, USA) placed inside a flow cabinet. The concentrations
of~selected elements (Ca, K, Mg, Na, P, S,~B, Mo, Cu, Fe, Mn, Zn) were
determined by inductively coupled plasma optical emission spectroscopy
ICP-OES (Thermo Jarrell-Ash, model 61E Polyscan, England).

\emph{Genomic analysis}

Leaf material from 2 A1\textsubscript{(c+)} individuals (F0-tolerant
pool), 2 T6\textsubscript{(c-)} individuals (F0-sensitive pool) and 19
F3 individuals from 7 tolerant families (2 T6xA1 and 5 A1xT6)
(F3-tolerant pool) were used for DNA extraction (Dataset SD5). DNA was
isolated with a Qiagen DNeasy kit and libraries were prepared using a
TruSeq DNA PCR-Free Sample Preparation Kit from Illumina. Whole-genome
sequencing was performed on Illumina HiSeq 2000 by Novogene (UK). Read
depth was measured after processing and alignment to the TAIR10
reference (below) at 15x coverage for each individual sample and 180x
for the pooled sample (19 \emph{F3} individuals).

Raw sequence data was processed as follows: (1) read trimming by default
parameters with Sickle {[}23{]}; (2) alignment to the TAIR10 reference
with bwa mem 0.7.12 (processing with samtools 1.3); (3) removal of
duplicate reads using Picard (MarkDuplicates); (4) application to read
group name correction to the bam files using Picard
(AddOrReplaceReadGroups); (5) realigning Indels using the GATK
`GenomeAnalysis' Toolkit {[}24{]}.

Bi-allelic SNPs were identified using `HaplotypeCaller' and genotyped
using `GenotypeGVCF' (both tools in GATK). Data were filtered using GATK
SelectVariants using these `GATK best practices' guidelines: QD
\textless{} 2.0 \textbar{}\textbar{} MQ \textless{} 40.00
\textbar{}\textbar{} FS \textgreater{} 60.0 \textbar{}\textbar{} SOR
\textgreater{} 4.0 \textbar{}\textbar{} MQRankSum \textless{} -8.0
\textbar{}\textbar{} ReadPosRankSum \textless{} -8.0 and a minimum
coverage of 10 per sample. Genotyping of the pooled sample
(19\emph{A1xT6} individuals) was performed with LoFreq {[}25{]}, as GATK
failed with ploidy levels required (ploidy-38 for 19 diploid samples).
LoFreq was run with default settings. VCFs were then merged with
BCFtools {[}26{]}. The allele frequency of all variants in the
F3-tolerant pool \emph{versus} the F0-sensitive pool was determined
using the Multipool program {[}27{]}. Several intervals with bi-polar
peaks were identified across the genome. Intervals where the allele
frequency between the two pools differed considerably and had a positive
LOD score were further explored for candidate genes.

The whole genome of 4 A1\textsubscript{(c+)} and 4
T6\textsubscript{(c-)} individuals was sequenced in a previous study
{[}28{]} and were merged with the 2 parental samples sequenced for the
BSA-seq experiment. To obtain a consensus sequence for each deme,
private SNPs of T6 and A1 (shared for the 6 samples) were selected using
GATK and Col-0 TAIR10 sequence as a reference. Once we obtained the
consensus sequence, we selected the sites that differ between demes
using GATK `concordance' command obtaining one vcf file with T16
AF\textless{}0.1 and A1 AF \textgreater{} 0.9 variants and a second vcf
file with T16 AF\textgreater{}0.9 and A1 AF\textless{} 0.1 variants. We
merged the files using VCFtools {[}29{]} and the amino acid changes
between A1 and T6 were obtained and quantified using SNPeff {[}30{]}.
Genes with 3 to 10 genetic modifiers were selected as candidates
(Dataset SD3).

\emph{RNA isolation and microarray scanning~~}

Root tissue from 15-day-old plants cultivated under hydroponic
conditions was recollected at 48 hours after initiating the
treatments.~At least 12 plants per line and treatment were pooled to
perform 3 biological replicates. Leaves were immersed in liquid
nitrogen, homogenized to a fine powder, and stored at -80o C. The total
RNA of about 100 mg of leaf powder for each biological replicate was
extracted using the Maxwell(r) plant RNA kit (Promega Corporation, WI,
USA) following the manufacturer's instructions. To ensure a high quality
and quantity RNA material, Experion RNA Analysis Kit (Biorad) was
performed. After quality control, total RNA was tag, hybridized and
cleaning using~ GeneChip(r) 3' IVT Express Kit. Microarray expression
analysis were performed by Affymetrix GeneChip\selectlanguage{english} Arabidopsis Genome ATH1
Array~Agilent 2100 Bioanalyzer.

\emph{Differential expression analysis:~} The function
\emph{gcrma}(package gcrma) was used to adjust background intensities in
Affymetrix array data, which include optical noise and non-specific
binding. The R package `limma' was used to test the expression data and
search for differentially expressed genes (DEGs) when pairs of
experimental groups were compared. A model matrix was defined by
specifying the groups that belonged to each sample. Then, limma's
functions \emph{lmFit} and\emph{contrasts.fit} , together with the model
matrix, were used to adjust the expression data to a linear model in
order to extract fold changes and confidence statistics associated to
the comparisons of interest. Function \emph{eBayes} was then used to
rank genes in order of evidence for differential expression.~ The R
library\emph{ath1121501.db} was used to include names and descriptions
of the genes associated to the microarray's features probed. The
resulting p-values were adjusted using the Benjamini and Hochberg's
approach {[}31{]} for controlling the false discovery rate (FDR). To
perform a list of DEGs, genes were filtered by~adj p-value (adjusted
p-value) \textless{}0.05 and LFC \textgreater{}1 and LFC \textless{}-1.

\emph{Gene ontology, KEGG pathway and functional protein association
networks analysis:} Gene Ontology (GO) enrichment analysis of DEGs was
implemented by AgriGO V2 {[}32{]} and FLAME {[}33{]}. Significant GO
terms (p-values less than 0.05) were classified into 3 categories:
biological function, molecular process, and cellular component. KEGG
pathway and STRING version 11.0 were used to understand high-level
function and gene interaction network from differential expressed genes
{[}34{]}.

\emph{Statistical analysis~of phenomics and ionomics}

Normal distribution of data was confirmed by Levene's test and
non-normal data were transformed before applying any parametric tests.
Mean-standardized values (1 \textless{} value \textgreater{}1) of
elemental contents of root material were used to represent the radar
plots and compare each accession profile. One-way or multivariate ANOVA
was used to test for significant differences (p value \textless{} 0.05)
between means of data with respect to physiological parameters, ionome,
fitness and gene expression. To test for correlations between two
variables, a bivariate fit was applied. To perform multiple comparisons
of group means, Tukey's HSD test was conducted. To perform multiple
comparisons to a control (WT), Dunnett's test was applied. Statistical
analyses and plots were performed using SAS Software JMP v.16.0 and are
available at Datasets SD2 and SD5.

\textbf{Results and discussion}

Natural variation of wild species provides a relevant complementary
resource to discover novel gene functions. It has been proven
that\emph{A. thaliana} demes from the north-east of Spain harbor high
genetic diversity {[}35{]} and they have a great potential for studying
traits related to adaptation (e.g. {[}36{]}) and allelic variants that
specifically interact with the environment (e.g. {[}37{]}).
Teres\emph{et al.} {[}5{]} selected two demes with contrasting
phenotypes in soil carbonate tolerance, A1\textsubscript{(c+)} and
T6\textsubscript{(c-),} and performed several studies to characterize
the physiological mechanisms activated by these demes under carbonate
and iron-deficiency stress.

\emph{Root characterization of A.thaliana demes under high pH and
bicarbonate stress}

The hydroponic cultivation of
A1\textsubscript{(c+),}T6\textsubscript{(c-)} and Col-0 plants under
control pH (pH 5.9), high pH (pH 8.3), or bicarbonate (bic, 10 mM
NaHCO\textsubscript{3}, pH 8.3) conditions confirmed the higher
tolerance of A1\textsubscript{(c+)} to both pH 8.3 and bic treatments.
Col-0 exhibited an intermediate response, while T6\textsubscript{(c-)}
was the most sensitive (Figure 1 \& {[}38{]}). Root length was hardly
affected by exposure to bic in A1\textsubscript{(c+)}, the deme evolved
on the calcareous soil. Contrastingly, root growth in
T6\textsubscript{(c-)} was strongly inhibited by bic, but less affected
by high pH (Figure 1A). In \emph{A. thaliana} Col-0, used as a
well-established reference genotype, root length was more inhibited by
bicarbonate than A1\textsubscript{(c+)}, but less than the sensitive
deme T6\textsubscript{(c-)}. In contrast to both natural demes, Col-0
did not show difference in root growth inhibition between high pH and
the bic treatment (Figure 1A).

The standardized root mineral content revealed clear differences in
response to high pH and bic among the three \emph{A. thaliana} variants
(Figure 1B). In deme A1\textsubscript{(c+)} both high pH and bic reduced
root Fe, Cu, and Zn concentrations in comparison to the control
treatment. However, bicarbonate stress enhanced the uptake of
macronutrients such as Ca and P in A1\textsubscript{(c+)} and Col-0
(Figure S1). In general, high pH and bic reduced Fe contents in both
organs in all plants. However, A1\textsubscript{(c+)} distinctively
maintained higher Fe levels under bicarbonate exposure and the Fe
translocation from roots to shoots was clearly enhanced in comparison to
the control treatment (Figure 1C). Contrastingly, the sensitive
T6\textsubscript{(c-)} had the lowest leaf Fe content despite having
higher root Fe levels (Figure S1). Because of this,
T6\textsubscript{(c-)} exhibited a clear inhibition by bicarbonate of Fe
translocation from roots to shoots, while high pH strongly enhanced the
Fe translocation in this sensitive deme (Figure 1C).

\emph{Transcriptional changes in response to bicarbonate stress}

Microarray analyses in roots were performed to characterize, at the gene
expression level, the mechanisms underlying the differential response to
bicarbonate of the two naturally selected accessions. Gene responses
were assessed by using samples after 48 h of treatment exposure, when
the plants still did not exhibit any foliar symptoms of the stress
treatments. We performed pairwise comparation to understand the
differences between A1\textsubscript{(c+)}-tolerant and
T6\textsubscript{(c-)} -sensitive lines in response to different
treatments (pH 8.3 \emph{vs} pH 5.9 and bic \emph{vs} control pH 5.9)
(Dataset SD3). A total of 4595 differentially expressed genes (DEGs)
were identified, considering accession and treatment (Figure 2). The
bicarbonate treatment caused a higher number of DEGs than the high pH
treatment, indicating that bicarbonate stress involves more complex
processes than simply the specific responses to alkaline pH (Figure 2A).
For that, we further focused our analysis on the pathways found in the
plants exposed to the bic treatment.

Previously, RNA-Seq of leaves from A1\textsubscript{(c+)} and
T6\textsubscript{(c-)} individuals revealed that bicarbonate exposure
quickly upregulated Fe-deficiency related genes in the sensitive
T6\textsubscript{(c-)} but not in the tolerant
A1\textsubscript{(c+)}{[}38{]}. In leaves, the highest number of DEGs
was observed in T6\textsubscript{(c-)} but in roots the tolerant
A1\textsubscript{(c+)}exhibited a huge response when exposed to bic
(Figure 2B). These differences in DEG number may denote contrasted deme
strategies toward bicarbonate stress. Gene Ontology (GO) terms indicated
that in roots of A1\textsubscript{(c+)} the main pathways activated
after 48 h of exposure were related to activation of biological
processes and metabolism located in the extracellular region (Figure 2C,
Dataset SD4). Also, enzymatic regulatory activity, transport and
transcription factors were altered. In T6\textsubscript{(c-)} only up
and down modifications in catalytic activity were found. The analysis of
the KEEG pathways revealed that the sensitive line activated only two
principal pathways involved in catalytic activity while in
A1\textsubscript{(c+)} more than 20 pathways were up or down regulated
(Figure 2D, Dataset SD4). This analysis indicates that the modifications
are being produced in the sugar, lipid and protein metabolism of the
tolerant line.

In rice, better tolerance to saline-alkaline conditions was due to
superior Fe acquisition and higher root to shoot Fe translocation
{[}39{]}. This was attributed to enhanced root development and higher
expression of genes involved in Fe acquisition and transport in the
tolerant variety in comparison to the sensitive one. Contrastingly,
among the DEGs in our tolerant \emph{A. thaliana} deme, no differential
enhancement of main Fe-acquisition related genes was observed. Moreover,
root expression of \emph{bHLH101} was more than 4 times higher in the
sensitive T6\textsubscript{(c-)} than in the tolerant
A1\textsubscript{(c+)} (Dataset S5). bHLH101 is a transcription factor
that controls plant Fe homeostasis in a FIT independent way {[}40{]}. It
is required for proper growth under Fe deficiency conditions. These
results further confirm that Fe deficiency is earlier perceived by the
sensitive T6\textsubscript{(c-).}

Focusing on A1\textsubscript{(c+)} \emph{vs}
T6\textsubscript{(c-)}comparison under the bic treatment and excluding
the genes also up or down regulated in control conditions, 208 DEGs were
identified in roots, 2503 DEGs in leaves, and 26 DEGs in both tissues
(Figure 2D, Dataset SD3). Two contrasted root strategies imply the
possibility of differences in signal perceiving and transduction. To
further explore these differential mechanisms, we performed
protein-protein interaction network functional enrichment analysis with
the STRING database on the 208 genes (Figure 3, Dataset SD4). These
genes can be classified in 5 main keywords according to FLAME:
catalytic, binding and ubiquitin activity; phosphorus metabolic process;
and oxidoreductase activity. We found two main cores of genes with high
network connectivity: one lead by genes involved in nucleus binding
(TOPII\textbf{,} CDC20.2, MAD2, AUR1, FZR3, CYCB1;3) and microtubule
movement regulation (AT5G60930, ATK5, AT3G20150, AT4G15830) that are
indispensable for normal plant development and fertility but have not
previously been associated with abiotic stress resistance; and another
controlled by the NIA1 and NIR1 genes. NIA1 (or GNR1) encodes the
cytosolic minor isoform of nitrate reductase and Nitrite Reductase 1
(NIR1) is involved in the second step of nitrate assimilation.
Contrastingly, \emph{NRT2} .3 was less expressed in
A1\textsubscript{(c+)} than in the sensitive T6\textsubscript{(c-)}.
NRT2.3 is a high affinity nitrate transporter with still poorly defined
functions but linked to the signaling of several phytohormones {[}41{]}.
The occurrence of inorganic carbon and nitrogen in karst and carbonated
soils affects the carbon/nitrogen metabolism of plant species. Under
bicarbonate stress, growth reduction is enhanced -due to the inhibition
of the photosynthesis and nitrogen metabolism- but water use efficiency
is promoted in tolerant plants {[}42{]}. This balance is vital for
plants to adapt to alkaline environments and the regulation of the
`nitrogen hub genes' found here might be important for C fixation and
bicarbonate tolerance.

Phenotypic variation for soil carbonate tolerance traits among demes

The natural habitat of A1\textsubscript{(c+)} and
T6\textsubscript{(c-)}differs mainly in soil pH (A1\textsubscript{(c+)}
: 7.4; T6\textsubscript{(c-)} : 6.5) and soil carbonate content
(A1\textsubscript{(c+)} : 12\%; T6\textsubscript{(c- )} : 0.8\%)
{[}5{]}. Natural carbonated soil from a location closer to
A1\textsubscript{(c+)} deme (pH = 7.9; CaCO\textsubscript{3} = 33\%) was
excavated for studying the tolerance capability of each deme and their
progeny. F1 progeny of both crosses (A1xT6 \& T6xA1) exhibited similar
growth and fitness as the tolerant parental
A1\textsubscript{(c+)}(Figure 4A\&C), while F2 progeny showed
segregation (Figure 4B\&C). In addition to symptoms of iron deficiency
and reduced growth, T6\textsubscript{(c-)} grown on alkaline soil
suffers from delayed flowering and infertility that substantially hamper
its fitness (Figure S2, Dataset SD5). Only the 35\% of F2 progeny plants
were able to flower and reproduce (Figure 4B). In consequence the
fitness of both crosses did not strictly fulfill the Mendelian phenotype
ratio (3:1) (Figure 4B\&C), suggesting that the `soil carbonate
tolerance phenotype' in our demes is partial-dominant or a polygenic
trait. Moreover, considering only the plants that were able to flower,
seed production was higher in the A1xT6 than in the T6xA1 offprings
(Figure 4C) pointing to a potential parental effect.

\emph{Identification of genes associated with soil carbonate tolerance
by BSA-Seq analysis}

During meiosis, recombination reasserts the complement of alleles
segregating in hybrid progeny {[}43{]}. Thus, the study of
F\textsubscript{2}~populations from contrasted phenotype crosses can
provide valuable information on genetic mechanism by association
{[}44{]}. With the purpose of performing a BSA-Seq study, plants from 65
families of the A1xT6 cross and from 28 families of the T6xA1 cross were
self-pollinated. F3 progeny plants were cultivated and screened under
the same soil and growth conditions to select the extreme pools (Dataset
SD5). We realized that the most sensitive families had very low
germination rates, or the seedlings died after a few days post
germination (NP plants), therefore we were not able to sequence the
sensitive pool. Despite this, considering the families with a higher
number of flowered individuals and with an elevated silique production,
7 families were selected as the tolerant pool (F3\textsubscript{TP})
(Figure 4D).

The allele frequency (AF) comparison of the F3\textsubscript{TP} with
the sensitive parental (T6\textsubscript{(c-)}) enabled the
identification of the top 0.5\% outlier SNPs from the empirical
distribution (Figure 4E). The 5617 SNPs of both tails were associated to
1119 genes, mainly located in chromosome 5. Within this chromosomal
group, the higher proportion of SNPs were assigned to the short arm
(Figure 4F, Dataset SD6). Salome et al. {[}43{]} found that crossovers
and segregation distortion in F2 populations from different \emph{A.
thaliana} accessions were more frequent in chromosomes 1 and 5. Here,
68\% of our phenotype-associated divergent outlier genes belong to chr5,
suggesting that this chromosome harbors key polymorphisms that
facilitate soil carbonate tolerance.

Additionally, we sequenced the whole genome of 6
A1\textsubscript{(c+)}and 6 T6\textsubscript{(c-)} individuals collected
from their natural habitat. It is well known that amino acid
substitutions in a protein can cause a drastic phenotypic effect
{[}45{]}. The AF comparison between A1\textsubscript{(c+)} and
T6\textsubscript{(c-)} picked out 977 genes with 3 to 10 non-synonymous
variants with high or moderate predicted effect (Dataset SD7). There are
many examples of genomic variants whose frequencies are correlated with
environmental variables and temporal changes consistent with natural
selection in \emph{A. thaliana} {[}35, 47{]}. A1\textsubscript{(c+)} is
locally adapted to carbonated soils, thus it was important to consider
the genomic differences between A1\textsubscript{(c+)} and
T6\textsubscript{(c-)} in order to narrow the BSA results down. The 69
genes in common between these two analyses were selected as the
``genomic candidate genes list'' (Figure 4G, Dataset SD8).

\emph{Association Analysis for prediction of candidate genes}

It is not feasible to determine the consequences of the missense
variants on gene expression, protein structure and function of an
elevated number of putative candidate genes. Alternatively, candidate
gene lists can be bounded conducting association analysis such as the
successfully applied to identify genes associated with plant
architecture in \emph{Brassica napus} {[}48{]}, water stress responses
in wheat {[}49{]}, or chilling stress tolerance in rice {[}50{]}. Here,
we performed an association analysis by combining the genomic and
transcriptomic data described above. The candidate gene lists obtained
from the BSA-Seq and the parental comparison analysis (Dataset SD6-8)
was compared with the root and leaf DEGs exclusive of the bic treatment
in A1\textsubscript{(c+)} \emph{vs} T6\textsubscript{(c-)} comparison
(Dataset SD3), resulting in 18 matching genes (Figure 5A\&B). The
expression pattern of these genes confirms that generally they are more
expressed in the roots of A1\textsubscript{(c+)} but also in the leaves
of the sensitive T6\textsubscript{(c-)} deme (Figure 5C\&D, Dataset
SD9). This suggests that roots of the tolerant A1\textsubscript{(c+)}are
reacting quickly to the adverse alkaline condition avoiding secondary
effects in the leaves. Contrastingly, in the sensitive
T6\textsubscript{(c-)} leaves are promptly stressed, especially due to
the inhibition of Fe translocation (Figure 1C), and leaf gene expression
adjusts to these alterations in order to alleviate further injury.

The SNP changes present in these 18 genes were explored and several
missense variants were identified in both demes. The high number of SNPs
present in AHH and PCH1 of A1\textsubscript{(c+),} and in INV-E of
T6\textsubscript{(c-)} stands out (Dataset SD9, Figure S4). Some of
these variants could be responsible for the observed gene expression
differences or be linked to structural changes that affect the
regulation of the gene. In order to evaluate the consequences of the
no-expression of these genes under carbonate stress, the knock-out
mutants of the 18 genes (T-DNA lines, Dataset SD1), together with the
wild type of Col-0 and A1\textsubscript{(c+)} and
T6\textsubscript{(c-)}, were cultivated on the same carbonated soil (LP)
and under the same control (pH 5.9) and bicarbonate hydroponic
conditions (10 mM NaHCO\textsubscript{3}, pH 8.3) used in our previous
studies.

Seeds of m03 (Sulphotransferase12) did not germinate in any condition
and therefore \emph{sot12} could not been tested. Instead,
\emph{gdh2}(m08) germinated normally under carbonated soil and we
observed an enhanced fitness of this mutant in this condition (Figure
S3A\&B). Both mutants \emph{csld2} (m15) and \emph{fbn6} (m16) coding,
respectively, for Cellulose Synthase Like D2 and Fibrillin6 died only
when treated with bicarbonate stress under hydroponic conditions (Figure
S3A\&B), suggesting hypersensitivity to this stress. CSLD2, is required
for normal root hair development {[}51{]}. However, the inhibition of
root hair development alone may not cause plant death. The failure to
survive is more likely related to the putative role of CSLD2 in
collaboration with CSLD3 and CSLD5 in hemicellulose biosynthesis
{[}52{]} affecting early plant development. FBN6 is also required for
normal plant development and \emph{fbn6} mutants have stunted root
growth. This has been related to alteration in sulphate reduction,
enhanced glutathione, and cadmium tolerance {[}53{]}. We have previously
shown that bicarbonate tolerance in \emph{A. thaliana} was associated
with the upregulation of leaf genes related to sulphate acquisition
{[}38{]}. So far FBN6 has mainly been studied in relation to
chloroplasts, while the role in root development is still unknown.
Clearly, further characterization of the role of FBN6 in plant
development under bicarbonate is required.

On the other hand, \emph{dao1} (m01), \emph{aah} (m07), and
\emph{jaz10}(m10) grew significantly better than the WT under bic stress
(Figure S3A) but only \emph{dao1} exhibited higher fitness when it was
cultivated on carbonated soil (Figure S3B). \emph{DAO1,} coding for a
2-oxoglutarate and Fe (II)-dependent dioxygenase, is a main player in
irreversible auxin degradation {[}54{]} and a higher production of
adventitious roots has been shown in the \emph{dao1-1} loss of function
mutants {[}55{]}. Downregulation of \emph{DAO1} could enhance root
development and Fe acquisition, being beneficial in high pH habitats.
Allantoate Amido Hydrolase (AHH) catalyzes the degradation of allantoin
yielding CO\textsubscript{2} and four molecules of
NH\textsubscript{4}{[}56{]}. High levels of allantoin and downregulation
of allantoin degradation has been found to enhance salt stress tolerance
in \emph{A. thaliana} {[}57{]}. Our findings regarding the bicarbonate
tolerance of the \emph{ahh} mutant (Figure S3A) and the lower root
expression of AHH in A1\textsubscript{(c+)} indicate (Figure 5C\&D) a
role for allantoin also in bicarbonate tolerance. In~\emph{A. thaliana}
, JAZ proteins act as repressors of JA signaling and plants
overexpressing \emph{AtJAZ10}are insensitive to JA {[}58, 59{]}. In our
case, \emph{jaz10} exhibited a better response to bicarbonate stress
than the WT and \emph{JAZ10} was less expressed in the tolerant
A1\textsubscript{(c+)} deme, suggesting that a downregulation of this
gene could trigger JA signals important for tolerant responses.

Trichome Birefringence Like 19 -\emph{tbl19} (m11) did not grow or
produced more siliques than Col-0 but were not affected by the carbonate
treatment, exhibiting higher tolerance (Figure S3A\&B). In both
demes\emph{TBL19} was downregulated in roots and upregulated in leaves
(Figure 5). Under bicarbonate exposure, the knock-out mutant had high
root biomass and the highest fitness as expressed in number of siliques
(Figure S3). TBL19, also known as AtXYBAT1 {[}63{]}, has
O-acetyltransferase activity and is probably involved in the acetylation
of cellohexaose, a hemicellulose in \emph{A. thaliana} cell walls. Root
apoplastic iron in \emph{A. thaliana} is mainly stored in the
operationally defined HC1 fraction of hemicellulose and this bound-Fe
can be mobilized by coumarin-type phenolic exudation {[}60{]}. Under
bicarbonate exposure, the tolerant A1\textsubscript{(c+)} has higher
exudation of these phenolics than the sensitive
T6\textsubscript{(c-)}{[}5{]}. The decrease of acetylation may enhance
Fe binding to cell wall hemicellulose and its accessibility to the
coumarin-type phenolics, so enhancing the root to shoot Fe translocation
specifically in A1\textsubscript{(c+)}.

Alkaline/Neutral Invertase (INV-E) encodes a chloroplast-targeted
protein, although it is expressed in root tissues in \emph{A.
thaliana}(TAIR-BAR). INV-E catalyzes the irreversible cleavage of
sucrose into glucose and fructose and some \emph{At} INV have been
linked to the control of root cell elongation mediated by sugars
{[}61{]}.\emph{inv-e} (m17) mutant was clearly affected by the
bicarbonate content both in hydroponics and in soil (Figure S3A\&B). The
SNP changes detected in A1\textsubscript{(c+)} and
T6\textsubscript{(c-)} were evaluated and we found that
T6\textsubscript{(c-)} INV-E gene contains a large number of variants
that could affect the expression of the gene (Figure S4, Dataset SD9),
supporting the hypothesis that a proper regulation of INV-E could play
an important role in bicarbonate tolerance.

In conclusion, differential bicarbonate tolerance in \emph{A.
thaliana}is associated with enhanced root to shoot translocation of
iron, proper distribution of C/N fixation, and better root development.
Association analysis of genomics and transcriptomics with natural demes
well- or mal-adapted to soil carbonate allowed the selection of 18
candidate genes for this differential tolerance. Evaluation of genomic
variation and the corresponding knock-out mutants suggested implications
of AHH, CSLD2, DAO1, GDH2, FBN6, INV-E, JAZ10 and TBL19. We propose
these genes as candidates for being further explored through
complementary analysis under bicarbonate stress conditions.

\textbf{Acknowledgments}

Special thanks to Sian Bray for her bioinformatic programming support.
We are grateful to Maria Jose Almira for her collaboration in laboratory
experiments and graphics; and to Rosa Padilla for processing the ionome
data. This research was funded by the Spanish Ministry of Science \&
Innovation (PID2019-104000RB-I00).

\textbf{References}

\begin{enumerate}
\tightlist
\item
  Li, N; Liu, H.; Sun, J.; Zheng, H.; Wang, J.; Yang, L.; Zhao, H.; Zou,
  D. Transcriptome analysis of two contrasting rice cultivars during
  alkaline stress. Sci. Rep. \textbf{2018} , 8, 9586.
  https://doi.org/10.1038/s41598-018-27940-x.
\item
  Filippi, P.; Jones, E.J.; Ginns, B.J.; Whelan, B.M.; Roth, G.W.;
  Bishop, T.F.A. Mapping the depth-to-soil pH constraint, and the
  relationship with cotton and grain yield at the within-field scale.
  Agronomy \textbf{2019} , 9, 251.
  https://doi.org/10.3390/agronomy9050251.
\item
  Loeppert, R.H.; Suarez, D.L. Carbonate and gypsum. In \emph{Methods of
  Soil Analysis. Part 3 Chemical Methods, 5.3} , Sparks, D.L.; Page,
  A.L.; Helmke, P.A.; Loeppert, R.H.; Soltanpour, P. N.; Tabatabai, M.
  A.; Johnston, C. T.; Sumner, M.E., Eds.; Publisher: Soil Science
  Society of America and American Society of Agronomy, Madison,
  USA,\textbf{1996} ; chapter 15, pp. 437-474, ISBN 978-089-118-825-4.
\item
  Al-Busaidi, A.; Cookson, P. Salinity -- pH relationships in calcareous
  soils of Oman. J. Sci Res. Agric. Marine Sci. \textbf{2003} , 8, 41 -
  46.
\item
  Teres, J.; Busoms, S.; Perez Martin, L.; Luis\selectlanguage{english}-Villarroya, A.; Flis,
  P.; Alvarez\selectlanguage{english}-Fernandez, A.; Tolra, R.; Salt, D.E.; Poschenrieder, C.
  Soil carbonate drives local adaptation in \emph{Arabidopsis thaliana}
  . Plant Cell Environ. \textbf{2019} , 42, 2384--98.
  https://doi.org/10.1111/pce.13567.
\item
  Hui, X.; Luo, L.; Wang, S.; Cao, H.; Huang, M.; Shi, M.; Malhi, S.S.;
  Wang, Z. Critical concentration of available soil phosphorus for grain
  yield and zinc nutrition of winter wheat in a zinc-deficient
  calcareous soil. Plant Soil \textbf{2019} , 444, 315-330,
  https://doi.org/10.1007/s11104-019-04273-w.
\item
  Dreyer, M.; Wichmann, M.; Rischen, M.; Gorlach, B.M.; Ehmke, A.;
  Pitmann, B.; Muhling, K.H. Ammonium-driven nitrification plays a key
  role in increasing Mn availability in calcareous soils. J. Plant Nutr.
  Soil Sci. \textbf{2020} , 183, 389-396,
  https://doi.org/10.1002/jpln.201900131.
\item
  Felle, H.H. pH: signal and messenger in plant cells. Plant
  Biol.\textbf{2001} , 3, 577-591. https://doi.org/10.1055/s-2001-19372
\item
  Lee, J.A.; Woolhouse, H.W. Root growth and dark fixation of carbon
  dioxide in calcicoles and calcifuges. New Phytol. \textbf{1969} , 68,
  247-255, https://doi.org/10.1111/j.1469-8137.1969.tb06437.x.
\item
  Sagervanshi, A.; Naeem, A.; Kaiser,H.; Pitann, B.; Muhling, K.H. Early
  growth reduction in \emph{Vicia faba} L. under alkali salt stress is
  mainly caused by excess bicarbonate and related to citrate and malate
  over accumulation. Environ. Exp. Bot. \textbf{2021} , 192, 104636,
  https://doi.org/10.1016/j.envexpbot.2021.104636.
\item
  Parvaiz, A.; Ozturk, M.; Sharma, S.; Gucel. S. Effect of sodium
  carbonate-induced salinity-alkalinity on some key osmoprotectants,
  protein profile, antioxidant enzymes, and lipid peroxidation in two
  mulberry (\emph{Morus alba} L.) cultivars. J. Plant
  Interact.\textbf{2014} , 9, 460--67,
  https://doi.org/10.1080/17429145.2013.855271.
\item
  Sarkar, B.; Hasanuzzaman, M.; Adak, M.K. Insights into the role of
  iron supplementation in conferring bicarbonate-mediated alkaline
  stress tolerance in maize. J. Soil Sci. Plant Nutr. \textbf{2022} ,
  22, 2719-2734, https://doi.org/10.1007/s42729-022-00839-3.
\item
  Zou, C.L.; Wang, Y.B.; Wang, B.; Liu, D.; Liu, L.; Li, C.F. Effects of
  alkali stress on dry matter accumulation, root morphology, ion
  balance, free polyamines, and organic acids of sugar beet. Acta
  Physiol. Plant. \textbf{2021} , 43, 13,
\item
  Cao, Y.; Song, H. and Zhang, L., 2022. New Insight into Plant
  Saline-Alkali Tolerance Mechanisms and Application to
  Breeding.~IJMS,~\textbf{2022} , 23, 16048.
  https://doi.org/10.3390/ijms232416048
\item
  Zhang, H..; Liu, X-L.; Zhang, R-X.; Yuan, H-Y.; Wang, M-M.; Yang,
  H-Y.; Ma, H-Y.; Liu, D.; Jiang, C.J.; Liang, Z-W. Root Damage under
  Alkaline Stress Is Associated with Reactive Oxygen Species
  Accumulation in Rice (\emph{Oryza sativa} L.). Front. Plant
  Sci.\textbf{2017} , 8, 1--12, https://doi.org/10.3389/fpls.2017.01580.
\item
  Waters, B.M.; Amundsen, K.; Graef, G. Gene expression profiling of
  iron deficiency chlorosis sensitive and tolerant soybean indicates key
  roles for phenylpropanoids under alkalinity stress. Front. Plant
  Sci.\textbf{2018} , 9, 10, https://doi.org/10.3389/fpls.2018.00010..
\item
  Duan, X.; Yu, Y.; Duanmu, H.; Chen, C.; Sun, X.; Cao, L.; Li, Q.;
  Ding, X.; Liu, B.; Zhu, Y. GsSLAH3, a \emph{Glycine soja} slow type
  anion channel homolog, positively modulates plant bicarbonate stress
  tolerance. Physiol. Plant. \textbf{2018} , 164, 145-162,
  https://doi.org/10.1111/ppl.12683.
\item
  Mehari TG, Xu Y, Umer MJ, Shiraku ML, Hou Y, Wang Y, Yu S, Zhang X,
  Wang K, Cai X, Zhou Z, Liu F. Multi-omics-based identification and
  functional characterization of GhA06G1257~proves its potential role in
  drought stress tolerance in~Gossypium \emph{hirsutum.} Front. Plant
  Sci. \textbf{2021} ,12, 746771,
  https:/doi.org/10.3389/fpls.2021.746771.
\item
  Li, Z.; Xu, Y. Bulk segregation analysis in the NGS era: a review of
  its teenage years.~Plant J. ~\textbf{2022} , 109, 1355-1374, .
\item
  Busoms, S. Teres, J.; Huang, X-Y.; Bomlies, K.; Danku, J.; Douglas,
  A.; Weigel, D.; Poschenrieder, C.; Salt, D.E. Salinity is an agent of
  divergent selection driving local adaptation of Arabidopsis to coastal
  habitats. Plant Physiol. \textbf{2015} , 168, 915-929.
  https://doi.org/10.1104/pp.15.00427
\item
  Scholl, R.L.; May, S.T.; Ware, D.H. Seed and molecular resources for
  Arabidopsis. Plant Physiol. \textbf{2000} , 124, 1477-1480.
  https://doi.org/10.1104/pp.124.4.1477.
\item
  Schneider, C.A. Rasband, W.S.; Eliceiri, K.W. NIH Image to ImageJ: 25
  years of image analysis. Nat. Methods 2012, 9, 671-5.
  https://doi.org/10.1038/nmeth.2089.
\item
  Joshi N.A.; Fass JN. Sickle: A sliding-window, adaptive, quality-based
  trimming tool for FastQ files (Version 1.33)
  {[}Software{]}.\textbf{2011.} Available at
  https://github.com/najoshi/sickle.
\item
  McKenna, A., et al.; DePristo, M.A., 2010. The Genome Analysis
  Toolkit: a MapReduce framework for analyzing next-generation DNA
  sequencing data.~Genome research,~\textbf{2010} , 20, 1297-1303.
  doi:~10.1101/gr.107524.110.
\item
  Wilm, A.; Aw, P.P.K.; Bertrand, D.; Yeo, G.H.T.; Ong, S.H.; Wong,
  C.H.; Khor, C.C.; Petric, R.; Hibberd, M.L.; Nagarajan, N. LoFreq: a
  sequence-quality aware, ultra-sensitive variant caller for uncovering
  cell-population heterogeneity from high-throughput sequencing
  datasets.~Nucleic acids research,~\textbf{2012} , 40, 11189-11201.
  https://doi.org/10.1093/nar/gks918.
\item
  Li, H. A statistical framework for SNP calling, mutation discovery,
  association mapping and population genetical parameter estimation from
  sequencing data.~Bioinformatics,~\textbf{2011} , 27, 2987-2993.
  https://doi.org/10.1093/bioinformatics/btr509
\item
  Edwards, M.D.; Gifford, D.K. High-resolution genetic mapping with
  pooled sequencing. BMC Bioinform. \textbf{2012} , 13 (Suppl 6), S8.
  https://doi.org/10.1186/1471-2105-13-S6-S8.
\item
  Busoms, S.; Paajanen, P.; Marburger, S.; Bray, S.; Huang, X.Y.;
  Poschenrieder, C.; Yant; Salt, D.E. Fluctuating selection on migrant
  adaptive sodium transporter alleles in coastal \emph{Arabidopsis
  thaliana} .~PNAS, \textbf{2018} ,~115, E12443-E12452.
  https://doi.org/10.1073/pnas.1816964115
\item
  Danecek, P., et al.; McVean, G. 1000 Genomes Project Analysis Group.
  2011.~\emph{The variant call format and VCFtools} .
  Bioinformatics,~\textbf{2011} , \emph{27} , 2156-2158.
  https://doi.org/10.1093/bioinformatics/btr330.
\item
  Cingolani, P. Variant Annotation and Functional Prediction: SnpEff.
  In: Ng, C., Piscuoglio, S. (eds) Variant Calling. Methods in Molecular
  Biology, \textbf{2022} , vol. 2493. Humana, New York, NY.
  https://doi.org/10.1007/978-1-0716-2293-3\_19.
\item
  Benjamini, Y.; Yekutieli, D. 2001, The control of the false discovery
  rate in multiple testing under dependency. Ann. Stat. \textbf{2001} ,
  29, 1165-1188. https://doi.org/10.1214/aos/1013699998.
\item
  Tian, T.; Liu, Y.; Yan, H.; You, Q.; Yi, X.; Du, Z.; Xu, W.; Su, Z.
  agriGOv2.0: a GO analysis for the agricultural community, 2017 update.
  Nucleic Acids Res. \textbf{2017} , 45, W122-W129.
  https://doi.org/10.1093/nar/gkx382.
\item
  Thanati, F.; Karatzas, E.; Baltoumas, F.; Stravopodis, D.J.;
  Eliopoulos, A.; Pavolpoulos, G. FLAME: a web tool for functional and
  literature enrichment analysis of multiple gene lists. Biol.
  (Basel)\textbf{2021} , 10:665.
  https://doi.org/10.3390/biology10070665.
\item
  Szklarczyk, D.; Gable, A.L.; Lyon, D.; Junge, A.; Wyder, S.;
  Huerta-Cepas, J.; Simonovic, M.; Doncheva, N.T.; Morris, H.; Bork, P.;
  Jensen, L.J.; von Mering, C. String v11: protein-protein association
  networks with increased coverage, supporting functional discovery in
  genome --wide experimental datasets. Nucleic Acids Res. \textbf{2019}
  , 47, D607-D613. https://doi.org/10.1093/nar/gky1131.
\item
  Castilla, A.R.; Mendez-Vigo, B.; Marcer, A.; Martinez-Minaya, J.;
  Conessa, D.; Pico, F.X.; Alonso-Blanco, C. Ecological, genetic and
  evolutionary drivers of regional genetic differentiation
  in\emph{Arabidopsis thaliana.} BMC Evol. Biol. \textbf{2020} , 20, 71,
  https://doi.org/10.1186/s12862-020-01635-2
\item
  Perez\selectlanguage{english}-Martin, L.; Busoms, S.; Almira, M.J.; Azagury, N.; Teres, J.;
  Tolra, R.; Poschenrieder, C.; Barcelo, J. Evolution of salt tolerance
  in \emph{Arabidopsis thaliana} on siliceous soils does not confer
  tolerance to saline calcareous soils. Plant Soil \textbf{2022} , 476,
  455--475, https://doi.org/10.1007/s11104-022-05439-9
\item
  Busoms, S.; Teres, J.; Yant, L.; Poschenrieder, C.; Salt, D.E.
  Adaptation to coastal soils through pleiotropic boosting of ion and
  stress hormone concentrations in wild \emph{Arabidopsis thaliana} .
  New Phytol. \textbf{2021} , 232, 208-220,
  https://doi.org/10.1111/nph.17569
\item
  Perez-Martin, L.; Busoms, S.; Tolra, R.; Poschenrieder, C.
  Transcriptomics reveals fast changes in salicylate and jasmonate
  signaling pathways in shoots of carbonate-tolerant \emph{Arabidopsis
  thaliana} under bicarbonate exposure. Int. J. Mol. Sci.,\textbf{2021}
  , 22, 1226.
\item
  Li, Q.; Yang, A.; Zhang, W-H. Efficient acquisition of iron confers
  greater tolerance to saline-alkaline stress in rice (\emph{Oryza
  sativa} L.). J. Exp. Bot. \textbf{2016} ,67, 6431--6444.
  https://doi.org/10.1093/jxb/erw407
\item
  Sivitz, A.B.: Hermand, V.; Curie, C.; Vert, G. Arabidopsis bHLH100 and
  bHLH101 control iron homeostasis via a FIT-independent pathway. PLoS
  One. \textbf{2012} , 7, e44843.
  https://doi.org/10.1371/journal.pone.0044843
\item
  Kiba, T.; Kudo, T.; Kojima, M.; Sakakibara, H. Hormonal control of
  nitrogen acquisition: roles of auxin, abscisic acid, and cytokinin. J.
  Exp. Bot. \textbf{2011} , 62, 1399-1409.
  https://doi.org/10.1093/jxb/erq410
\item
  Xia, A.; Wu, Y. Joint interactions of carbon and nitrogen metabolism
  dominated by bicarbonate and nitrogen in \emph{Orychophragmus
  violaceus} and \emph{Brassica napus} under simulated karst
  habitats.~BMC Plant Biol. \textbf{2022} , 22, 1-13.
  https://doi.org/10.1186/s12870-022-03646-1
\item
  Salome, P.A.; Bomblies, K.; Fitz, J.; Laitinen, R.A.E.; Warthmann, N.;
  Yant, L.; Weigel, D. The recombination landscape in \emph{Arabidopsis
  thaliana} F2 populations.~Heredity,~\textbf{2012} , 108, 447-455.
  https://doi.org/10.1038/hdy.2011.95
\item
  Benowicz, A.; Stoehr, M.; Hamann, A.; Yanchuk, A.D. Estimation of the
  F2 generation segregation variance and relationships among growth,
  frost damage, and bud break in coastal Douglas-fir (Pseudotsuga
  menziesii (Mirb.) Franco) wide-crosses.~Annals of Forest
  Science,~\textbf{2020} , 77, 1-13.
  https://doi.org/10.1007/s13595-020-0925-9
\item
  Bartlett, M.E.; Whipple, C.J. Protein change in plant evolution:
  tracing one thread connecting molecular and phenotypic
  diversity.~Front. Plant Sci.,~\textbf{2013} , 4, 382,
  https://doi.org/10.3389/fpls.2013.00382
\item
  Castilla, A.R.; Mendez-Vigo, B.; Marcer, A.; Martinez-Minaya, J.;
  Conessa, D.; Pico, F.X.; Alonso-Blanco, C. Ecological, genetic and
  evolutionary drivers of regional genetic differentiation
  in\emph{Arabidopsis thaliana.} BMC Evol. Biol. \textbf{2020} , 20, 71,
  https://doi.org/10.1186/s12862-020-01635-2
\item
  Monroe, J.G., Srikant, T., Carbonell-Bejerano, P.~\emph{et al.}
  ~Mutation bias reflects natural selection in~\emph{Arabidopsis
  thaliana} .~Nature, \textbf{2022} ,~602, 101--105.
  https://doi.org/10.1038/s41586-021-04269-6
\item
  Ye, S.; Yan, L.; Ma, X.; Chen, Y.; Wu, L.; Ma, T.; Zhao, L.; Yi, B.;
  Ma, C.; Tu, J.; Shen, J. Combined BSA-Seq based mapping and RNA-Seq
  profiling reveal candidate genes associated with plant architecture
  in\emph{Brassica napus} .~Int. J. Mol. Sci. ~\textbf{2022} , 23, 2472.
  https://doi.org/10.3390/ijms23052472
\item
  Derakhshani, B.; Ayalew, H.; Mishina, K.; Tanaka, T.; Kawahara, Y.;
  Jafary, H.; Oono, Y. Comparative analysis of root transcriptome
  reveals candidate genes and expression divergence of homoeologous
  genes in response to water stress in wheat.~Plants,\textbf{2020}
  ,~9596. https://doi.org/10.3390/plants9050596
\item
  Guo, Z.; Cai, L.; Chen, Z.; Wang, R.; Zhang, L.; Guan, S.; Zhang, S.;
  Ma, W.; Liu, C.; Pan, G. Identification of candidate genes controlling
  chilling tolerance of rice in the cold region at the booting stage by
  BSA-Seq and RNA-Seq. ~R. Soc. Open Sci.,~\textbf{2020} , 7.201081.
  http://dx.doi.org/10.1098/rsos.201081
\item
  Bernal, A.J.; Yoo, C-M.; Mutwil, M.; Jensen, J. K.; Hou, G.; Blaukopf,
  C.; Sorensen, I.; Blancaflor, E. B.; Scheller, H. V.; Willats, W.G.T.
  Functional analysis of the cellulose synthase-like genes CSLD1, CSLD2,
  and CSLD4 in tip-growing Arabidopsis cells. Plant
  Physiol.\textbf{2008} , 148, 1238--1253.
  https://doi.org/10.1104/pp.108.121939
\item
  Yin, L.; Verhertbruggen, Y.; Oikawa, A.; Manisseri, C.; Knierim, B.;
  Prak, L.; Kruger Jensen, J.; Knox, J.P.; Auer, M.; Willats, W.G.T.;
  Scheller, H.V. The cooperative activities of CSLD2, CSLD3, and CSLD5
  are required for normal Arabidopsis development. Mol.
  Plant\textbf{2011} , 4, 1024-1037. https://doi.org/10.1093/mp/ssr02
\item
  Lee, K.; Lehmann, M.; Paul, M.V; Wang, L.; Luckner, M.; Wanner, G.;
  Geigenberger, P.; Leister, D.; Kleine, T. Lack of FIBRILLIN6
  in\emph{Arabidopsis thaliana} affects light acclimation and sulfate
  metabolism. New Phytol. \textbf{2020} , 225, 1715-1731.
  https://doi.org/10.1111/nph.16246
\item
  Hayashi, K.; Arai, K.; Tanaka, Y.; Hira, H.; Guo, R.; Hu, Y.; Ge, C.;
  Zhao, Y.; Kasahara, H.; Fukui, K. The main oxidative inactivation
  pathway of the plant hormone auxin. Nat. Commun. \textbf{2021} , 12,
  6752, https://doi.org/10.1038/s41467-021-27020-1
\item
  Lakehal, A.; Dob, A.; Novak, O.; Bellini, C. A DAO1-mediated circuit
  controls auxin and jasmonate crosstalk robustness during adventitious
  root initiation in Arabidopsis. Int. J. Mol. Sci. \textbf{2019} , 20,
  4428, https://doi.org/10.3390/ijms20184428
\item
  Werner, A.K.; Sparkes, I.A.; Witte, C-P. Identification, biochemical
  characterization, and subcellular localization of allantoate
  amidohydrolasres from Arbidopsis and soybean. Plant
  Physiol.\textbf{2008} , 146, 418-430
  https://doi.org/10.1104/pp.107.110809
\item
  Lescano, C.I.; Martini, C.; Gonzalez, C.A.; Desimone, Allantoin
  accumulation by allantoinase downregulation and transport by ureide
  permease 5 confers salt stress tolerance to Arabidopsis plants. Plant
  Mol. Biol. \textbf{2016} , 91, 581-595,
  https://doi.org/10.1007/s11103-016-0490-7
\item
  Chung, H.S.; Howe, G.A. A critical role for the TIFY motif in
  repression of jasmonate signaling by a stabilized splice variant of
  the JASMONATE ZIM-domain protein JAZ10 in Arabidopsis. Plant
  Cell,\textbf{2009} , 21, 131-145.
  https://doi.org/10.1105/tpc.108.064097
\item
  Sun, Q.; Wang, G.; Zhang, X.; Zhang, X.; Qiao, P.; Long, L.; Yuan, Y.
  and Cai, Y., 2017. Genome-wide identification of the TIFY gene family
  in three cultivated Gossypium species and the expression of JAZ
  genes.~ Sci. Rep. ~\textbf{2017} , 7, 1-9.
  https://doi.org/10.1038/srep42418
\item
  Lei, G.J.; Zhu, X.F.; Wang, Z.W.; Dong, F.; Dong, N.Y.; Zheng, S.J.
  Abscisic acid alleviates iron deficiency by promoting root iron
  reutilization and transport from root to shoot in Arabidopsis. Plan
  Cell Environ. \textbf{2014} , 37, 852-863.
  https://doi.org/10.1111/pce.12203
\item
  Ruan, Y.L.; Jin, Y.; Yang, Y.J.; Li, G.J.; Boyer, J.S. Sugar input,
  metabolism, and signaling mediated by invertase: roles in development,
  yield potential, and response to drought and heat.~Mol.
  Plant,\textbf{2010} ,~3, 942-955. https://doi.org/10.1093/mp/ssq044
\end{enumerate}

\textbf{Tables and Figure legends}

\textbf{Table 1.} Physical and chemical properties of the natural
carbonated soil excavated from Les Planes d'Hostoles (LP, 42deg 03\selectlanguage{english}' 45.1"
\selectlanguage{english}N; 2\selectlanguage{ngerman}° 32\selectlanguage{english}' 46.6" \selectlanguage{english}E). O.M.= organic matter; mineral nutrient
concentrations are means \selectlanguage{ngerman}± SD (SD= standard deviation).

\textbf{Figure 1. Root growth and ionome. (A)} Root length
(cm),\textbf{(B)} root ionome profile (mean-standardized values of 12
elements), and \textbf{(C)} root Fe concentration (\selectlanguage{greek}μ\selectlanguage{english}g/g) and Fe
translocation ratio (leaf {[}Fe{]} / root {[}Fe{]}) of A1(c+), T6(c-)
and Col-0 \emph{A. thaliana} plants cultivated in hydroponics under
control (pH 5.7), high pH (pH 8.3) and bic (10 mM NaHCO3, pH 8.3)
conditions for 10 days. Numbers indicate significant differences between
treatments per accession and letters indicate significant differences
between accessions under the same treatment (t-test, p \textless{}
0.05).

\textbf{Figure 2} . \textbf{Transcriptomics. (A)} Total number and
heatmap profile of differentially expressed genes (DEGs) from pairwise
comparisons bic \emph{vs} control (dark blue) and high pH
\emph{vs}control (light blue) in A1\textsubscript{(c+)} and
T6\textsubscript{(c-)} demes. \textbf{(B)} Total number of pairwise
comparisons bic \emph{vs} control in roots and leaves of
A1\textsubscript{(c+)} and T6\textsubscript{(c-)} demes. Bubble plots
indicating \textbf{(C)} significant Gene Ontology (GO) analysis
and\textbf{(D)} heatmap of KEGG pathway analysis of
A1\textsubscript{(c+)}and T6\textsubscript{(c-)} root DEGs in bic
\emph{vs} control comparison. Arrows indicate up or downregulated genes.
\textbf{(E)} Venn diagrams of pairwise comparisons
A1\textsubscript{(c+)} bic \emph{vs}T6\textsubscript{(c-)} bic and
A1\textsubscript{(c+)} control \emph{vs}T6\textsubscript{(c-)} control
in roots and leaves. Selected DEGs are highlighted in brown for root and
green for leaves. DEGs, GO and KEGG terms were filtered at log fold
change (LFC \textgreater{} 1, LFC \textless{} -1), and adjusted p-value
\textless{} 0.05.

\textbf{Figure 3. Protein--protein interaction network of bicarbonate
stress} . Gene protein interaction network of the 208 root DEGs up or
down regulated in the tolerant A1\textsubscript{(c+)} under the bic
treatment. Each sphere corresponds to one gene, nodes represent protein
interactions (source:~STRING.org) and colors represent GO terms (source:
bib.fleming.gr).

\textbf{Figure 4.} \emph{Arabidopsis thaliana} \textbf{soil carbonate
tolerance inheritance and BSA-Seq} . \textbf{(A)} Rosette diameter
(mm);\textbf{(B)} percentage of flowered plants and plants with more or
less than 10 siliques; and \textbf{(C)} fitness (silique number) of
A1\textsubscript{(c+)}, T6\textsubscript{(c-)} parents, F1 and F2
progeny from reciprocal crosses cultivated in natural carbonated
soil.\textbf{(D)} Percentage of plants that did not prosper (NP, purple
bars), plants that did not flower (NF, grey bars), and plants that
produce more than 10 siliques (orange bars) from A1xT6 F3 and T6xA1 F3
progeny cultivated in natural carbonated soil. Selected families for BSA
analysis are indicated in bold. \textbf{(E)} Allele frequency
distribution of SNPs detected between the sensitive parent
(T6\textsubscript{(c-)}) and the F3 tolerant pool.
\textbf{(F)}Chromosome distribution of genes associated with the SNPs
from BSA left and right tails. \textbf{(G)} Summary table of the SNPs
and genes detected and shared between parent's comparison
(A1\textsubscript{(c+)}\emph{vs} T6\textsubscript{(c-)}) and between
T6\textsubscript{(c-)} and the F3 tolerant pool comparison. Letters
indicate significant differences (Tukey's HSD, p \textless{} 0.05).

\textbf{Figure 5. Combined analysis for bicarbonate tolerance candidate
genes detection. (A)} Number and \textbf{(B)} description of putative
candidate genes obtained from the genomic and transcriptomic association
analysis of A1\textsubscript{(c+)} and T6\textsubscript{(c-)} demes.
Relative fold change (bic \emph{vs} control) of the 18 candidate genes
in the \textbf{(C)} leaves and \textbf{(D)} roots of
A1\textsubscript{(c+)} (orange bars) and T6\textsubscript{(c-)} (purple
bars) plants submitted to bic stress (10 mM NaHCO\textsubscript{3}, pH
8.3) for 48 hours.

\textbf{Supporting Information}

\textbf{Figure S1. (A)} Root and leaf ion concentrations (\selectlanguage{greek}µ\selectlanguage{english}g/g)
and\textbf{(B)} ion translocation ratios (leaf {[}X{]} / root {[}X{]})
of A1\textsubscript{(c+)}, T6\textsubscript{(c-)} and Col-0 \emph{A.
thaliana} plants cultivated in hydroponics under control (pH 5.7), high
pH (pH 8.3) and bic (10 mM NaHCO\textsubscript{3}, pH 8.3) conditions
for 10 days. Letters indicate significant differences between treatments
for each accession (t-test, p \textless{} 0.05).

\textbf{Figure S2. (A)} Growth (rosette diameter, mm) and fitness
(silique number) correlation and stats. \textbf{(B)} Bolting time and
fitness correlation and stats.

\textbf{Figure S3. (A)} Aerial and root biomass (dry weight, g) of
A1\textsubscript{(c+)}, T6\textsubscript{(c-),}Col-0\textsubscript{(WT)}
and 18 T-DNA knockout mutants of the candidate genes grown
hydroponically under control (1/2 Hoagland Solution (HS), pH 5.9) and
bic stress (1/2 HS + 10 mM NaHCO\textsubscript{3}, pH 8.3) for 10 days.
NP=no plant (no germination in phytoagar); D=dead plants.\textbf{(B)}
Fitness (silique number) of A1\textsubscript{(c+)},
T6\textsubscript{(c-),} Col-0\textsubscript{(WT)} and 18 T-DNA knockout
mutants cultivated in natural carbonated soil. NP=no plant (no
germination in carbonated soil). Letters indicate significant
differences between treatment (p \textless{} 0.05, t-test) and numbers
or asterisks indicate significant differences respect the wild-type
Col-0 (p \textless{} 0.05, Dunnet's test).

\textbf{Figure S4.} SNP changes detected in the \textbf{(A)} AHH
(AT4G20070), \textbf{(B)} PCH1 (AT2G16365), and \textbf{(C)} INV-E
(AT5G22510) genes of A1\textsubscript{(c+)} and T6\textsubscript{(c-)}.
Missense variants are marked in red.

\textbf{Dataset S1. Plant materials.} \emph{A. thaliana} natural demes
information, T-DNA lines of 18 candidate genes and primers used for qPCR
validation.

\textbf{Dataset S2} . Phenomic and ionomic statistics from hydroponic
and carbonated soil cultivations.

\textbf{Dataset S3} . DEGs from A1\textsubscript{(c+)} and
T6\textsubscript{(c-)} leaf and root transcriptomic analysis.

\textbf{Dataset S4.} GO terms, KEGG pathways, and STRING output from
root transcriptomics.

\textbf{Dataset S5.} Fitness and screening of A1xT6 progeny under
natural carbonated soil.

\textbf{Dataset S6.} BSA-Seq outliers (1119 genes) and GO terms of
outliers from Chromosome 5.

\textbf{Dataset S7.} A1\textsubscript{(c+)} vs
T6\textsubscript{(c-)}allele frequency comparison (977 genes)\textbf{.}

\textbf{Dataset S8.} Input gene lists for genomic and transcriptomic
association analysis\textbf{.}

\textbf{Dataset S9.} Gene expression and SNP changes of the 18 putative
candidate genes.\selectlanguage{english}
\begin{figure}[H]
\begin{center}
\includegraphics[width=0.70\columnwidth]{figures/DI-Figure-1-foot/DI-Figure-1-foot}

\end{center}
\end{figure}\selectlanguage{english}
\begin{figure}[H]
\begin{center}
\includegraphics[width=0.70\columnwidth]{figures/DI-Figure-2-foot/DI-Figure-2-foot}

\end{center}
\end{figure}\selectlanguage{english}
\begin{figure}[H]
\begin{center}
\includegraphics[width=0.70\columnwidth]{figures/DI-Figure-3-foot/DI-Figure-3-foot}

\end{center}
\end{figure}\selectlanguage{english}
\begin{figure}[H]
\begin{center}
\includegraphics[width=0.70\columnwidth]{figures/DI-Figure-4-foot/DI-Figure-4-foot}

\end{center}
\end{figure}\selectlanguage{english}
\begin{figure}[H]
\begin{center}
\includegraphics[width=0.70\columnwidth]{figures/DI-Figure-5-foot/DI-Figure-5-foot}

\end{center}
\end{figure}\selectlanguage{english}
\begin{figure}[H]
\begin{center}
\includegraphics[width=0.70\columnwidth]{figures/DI-Table-1/DI-Table-1}

\end{center}
\end{figure}

\selectlanguage{english}
\FloatBarrier
\end{document}

