Didier F Pisani
email: didier.pisani@univ-cotedazur.fr

Gérard Ailhaud

Involvement of Polyunsaturated fatty acids in the Control of Energy Storage and Expenditure

Keywords: Thermogenesis, fever, shivering, PUFA, oxylipins

Within heterotrophic organisms, fat, sugar and protein are required to bring energy. In humans, energy homeostasis (i.e. the balance between energy intake and energy expenditure) is tightly regulated. Obesity, characterized by an excess of white fat mass, is a consequence of dysregulation of this balance in which decreased energy expenditure plays an important role. Among dietary components, fat represents approximately 30% of daily intake. Recent studies have shown that, besides its caloric input, fatty acid composition of fat represents an important qualitative issue. The impact of ω6 and ω3 polyunsaturated fatty acids, on weight gain and development has been extensively studied. Interestingly, the role played by these polyunsaturated fatty acids in energy expenditure has been also characterized and will be discussed in relation to the various organs involved, in particular brown adipose tissue.

Obesity, adipose tissue and nutrition

Obesity results from an imbalance between energy input and expenditure. This energy homeostasis is maintained by the activity of metabolic tissues, including white and brown adipose tissue. As adipose tissue is the main peripheral target organ handling fatty acids, any change in lipid intake involves modifications in its development and functions. Importantly, major dietary quantitative and qualitative changes have taken place in the last decades in the Western industrialized world, particularly the overconsumption of ω6 polyunsaturated fatty acids (PUFAs) and the underconsumption of ω3 PUFAs by both humans and domesticated animals, which have been associated to an increase in white fat mass leading to obesity. Over the last years, in addition to the well-known relationship between PUFAs and energy storage, new data have emerged linking PUFA to the different mechanisms of energy expenditure.

In this review we listed the various effects of ω6 and ω3 PUFA and long-chain PUFAs to explain that inadequate intake of PUFA contributes to energy homeostasis disequilibrium and may represent a risk factor to the occurrence of obesity and metabolic disorders.

Energy expenditure and thermogenesis

Energy expenditure involves resting and adaptative energy expenditure (Fig. 1).

Resting energy expenditure or resting metabolic rate includes physical activity (excluding exercise), diet-induced thermogenesis and all energy expenses required to sustain life functions (named basal metabolic rate and corresponding to breathing, blood circulation, brain and nerve functions). Adaptative energy expenditure is known as a response to specific situations, such as fever with shivering in response to infection and non-shivering thermogenesis in response to cold. Both processes, can be modulated/affected by exogenous stimuli which include dietary fatty acids [START_REF] Himms-Hagen | Does brown adipose tissue (BAT) have a role in the physiology or treatment of human obesity?[END_REF].

Non-shivering thermogenesis

Non-shivering thermogenesis is activated in response to a cold environment and corresponds to an increase in metabolic heat production above basal metabolism rate [START_REF] Cannon | Brown adipose tissue: function and physiological significance[END_REF]. Distinct from shivering thermogenesis, it is not exclusively associated with muscle activity. It occurs mainly through metabolism of brown fat and less to that of skeletal muscle, liver, brain, and white fat (Cannon andNedergaard 2004, Ouellet, Labbe et al. 2012).

In contrast to white adipose tissue (WAT) involved in energy storage and release, brown adipose tissue (BAT) is endowed with this thermogenic activity to regulate body temperature by dissipating energy. This process is due to the unique occurrence of the Uncoupling Protein 1 (UCP1) localized in BAT mitochondria and is both activated and induced in response to cold via β-adrenergic stimulation [START_REF] Cannon | Brown adipose tissue: function and physiological significance[END_REF]. The energy-dissipating properties of UCP1 lead to an increased oxidation of fatty acids and carbohydrates and are important for body weight regulation [START_REF] Blondin | Dietary fatty acid metabolism of brown adipose tissue in coldacclimated men[END_REF]. BAT displays various anatomic localizations, in humans where the interscapular localization found only in newborns is replaced by cervical and thoracic localization in adults [START_REF] Zhang | An Adipose Tissue Atlas: An Image-Guided Identification of Human-like BAT and Beige Depots in Rodents[END_REF]. In other mammals, the interscapular localization of BAT is maintained throughout adult life in addition to other sites, but its mass and activity decreases with aging [START_REF] Graja | Mechanisms of aging-related impairment of brown adipocyte development and function[END_REF].

Importantly, an additional population of UCP1-positive adipocytes is present in WAT and is termed brite for "brown in white" or beige adipocytes [START_REF] Barbatelli | The emergence of cold-induced brown adipocytes in mouse white fat depots is determined predominantly by white to brown adipocyte transdifferentiation[END_REF]). These brown-like adipocytes appear in response to cold exposure or to high-fat diets. These cells which stem from progenitors or by direct conversion of mature white adipocytes have been recently found in adult humans [START_REF] Rosenwald | The origin and definition of brite versus white and classical brown adipocytes[END_REF]. Induction of brown and brite adipocyte activity appears as a novel strategy to combat obesity by enhancing body energy expenditure, i.e.by increasing oxidation of fatty acids within these cells and therefore limiting fat storage in peripheral tissues.

Diet-induced thermogenesis

In addition to the well-known function of BAT in cold acclimation, several other functions of this tissue have been described. The most intriguing response is its postprandial activation in response to various diets, i.e. diet-induced thermogenesis which takes place concomitantly to the insulin peak which favors energy storage as triglycerides and glycogen [START_REF] Hibi | Brown adipose tissue is involved in diet-induced thermogenesis and whole-body fat utilization in healthy humans[END_REF].

Physical activity

"Physical activity" includes occupational, sports, conditioning, household, or other activities whereas "Exercise" is a subset of planned and structured physical activity. Importantly, physical activity is mainly performed by rapid skeletal muscles known as glycolytic and thus uses preferentially carbohydrates. Thus, it is likely that the fat composition of diet should have lower impact, if any, on this kind of energy expenditure [START_REF] Westerterp | Energy expenditure and physical activity in subjects consuming full-or reduced-fat products as part of their normal diet[END_REF].

Shivering thermogenesis

Exposure to cold induces the activation of involuntary muscle contractions leading to shivering thermogenesis. Skeletal muscles are thus involved in cold resistance and are the main contributor of heat in adult humans compared to other mammals in which BAT thermogenesis prevails [START_REF] Blondin | Dietary fatty acid metabolism of brown adipose tissue in coldacclimated men[END_REF]. In humans, both types of thermogenesis seem more exclusive than additive as shivering activity is inversely correlated to the occurrence of BAT [START_REF] Ouellet | Brown adipose tissue oxidative metabolism contributes to energy expenditure during acute cold exposure in humans[END_REF].

Shivering thermogenesis uses a variety of metabolic fuels, as oxidation of carbohydrates, lipids and proteins. When one fuel source is depleted or reduced, the others compensate in order to maintain shivering intensity and core temperature [START_REF] Weber | Fuel selection in shivering humans[END_REF].

Fever

Fever, or pyrexia, is an evolutionarily conserved host thermogenic response to microbial infection. It is defined as an elevated core temperature and is a reflect of an enhanced metabolism characterized by an increase in energy expenditure. Fever is mediated by pyrogenic mediators, such as cytokines (interleukin-1 and -6, IL-1 and IL-6) [START_REF] Kozak | IL-6 and IL-1 beta in fever. Studies using cytokine-deficient (knockout) mice[END_REF], adipokines (leptin) [START_REF] Luheshi | Leptin actions on food intake and body temperature are mediated by IL-1[END_REF] or eicosanoids derived from polyunsaturated fatty acids (prostaglandin E2, PGE2) [START_REF] Engblom | Microsomal prostaglandin E synthase-1 is the central switch during immune-induced pyresis[END_REF] which trigger thermogenic recruitment of peripheral tissues [START_REF] Evans | Fever and the thermal regulation of immunity: the immune system feels the heat[END_REF]. This response is counteracted by anti-pyrogenic agents (corticoids, arginin-vasopressin and IL-1 receptor antagonist, IL-1RA) which prevent an harmful increased body temperature [START_REF] Tatro | Endogenous antipyretics[END_REF]. In this respect, the involvement of BAT in fever remains disputed as, on one hand, inflammatory signals allow secretion of norepinephrine by central nervous system, which is considered as a major activator of BAT but, on the other hand, BAT seems excluded from the process which triggers fever in response to endotoxemia [START_REF] Riley | The complementary and divergent roles of uncoupling proteins 1 and 3 in thermoregulation[END_REF] .

Energy expenditure in response to PUFAs

Polyunsaturated fatty acids (PUFAs) are defined as a class of fatty acids displaying two or more double bonds. The two major PUFAs are ω6 linoleic acid (LA, 18:2n-6), the precursor of dihomo-γ-linolenic acid (DGLA, andarachidonic acid (ARA, 20:4n-6), andω3 α-linolenic acid (LNA, 18:3n-3), the precursor of ω3 eicosapentaenoic (EPA, 20:5n-3) anddocosahexaenoic (DHA, 22:6n-3 In animals, PUFAs are considered as essentials for health owing to their low efficiency (for DHA and to a lower extent for ARA and EPA) or to their inability (for LA, LNA) to synthesize them. Indeed, animals lack the Δ-12-and Δ-15-desaturases able to introduce double bonds in C18 long chain fatty acids. LA can be found at high proportions in most consumed oils [sunflower (75%), safflower (79%), evening primrose seed (72%), corn (57%), peanut (31%) and brassica (19%)]. ALA is found in linseed (also named flaxseed) and various oils [linseed (50%), black currant seed (14%) and for a lower extent in rapeseed (9%)]. EPA and DHA are present mainly in algae and marine bacteria. PUFAs can be also ingested via meat and fish, due to their high consumption of ω6 and ω3 PUFAs. Indeed, ARA is found in high amounts in bovine, porcine, ovine and poultry meats as well as in fish flesh, milk and eggs. Intake of EPA and DHA is less important than ARA intake as only found in wild fishes 'salmon, mackerel and tuna). In order to increase EPA and DHA intake in humans two strategies are developed, using dietary supplements (fish oil, seed oil) or supplementing livestock diet (especially cow, pork, chicken and salmon) with ω3-enriched food like linseed and algae [START_REF] Weill | Effects of introducing linseed in livestock diet on blood fatty acid composition of consumers of animal products[END_REF][START_REF] Bourre | Where to find omega-3 fatty acids and how feeding animals with diet enriched in omega-3 fatty acids to increase nutritional value of derived products for human: what is actually useful ?[END_REF].

The critical ω6 and ω3 PUFAs intakes in obesity and chronic diseases

Both ω6 and ω3 PUFAs are transported into the bloodstream and are incorporated as phospholipids in plasma membranes or as triglycerides within adipocyte droplets. PUFAs are released by lipases and then oxidized. ω6 and ω3

PUFAs are known to compete at different enzymatic steps that modulate the availability of their respective metabolic products (Fig. 2). LA and LNA are substrates of the same Δ-5-and Δ-6-desaturases and elongases to generate the LC-PUFAs ARA, DGLA, EPA and DHA. They are also competitors as substrates for the synthesis of oxygenated derivatives or for β-oxidation, ω3 PUFAs being more rapidly oxidized compared to ω6

PUFAs and to mono-unsaturated fatty acids [START_REF] Cunnane | Problems with essential fatty acids: time for a new paradigm?[END_REF]. Due to these different competition levels, the intake amount of ω6 and ω3 PUFAs is crucial to determine the kind of mediators available and the downstream biological processes which will be ultimately induced [START_REF] Bibus | Balancing proportions of competing omega-3 and omega-6 highly unsaturated fatty acids (HUFA) in tissue lipids[END_REF]. For a long-time, dietary recommendations were mainly addressed to pregnant women and based on an optimal LA/LNA ratio for a healthy development of fetuses and children [START_REF] Simopoulos | The importance of a balanced omega-6 to omega-3 ratio in the prevention and management of obesity[END_REF].

Numerous studies carried out in rodent models and in human cohorts have shown a tight link between increased ω6 PUFA intake, generally associated with a high LA/LNA ratio, and enhanced fat mass. In mice, high fat diet rich in ω6 PUFA favored fat accumulation, in contrast to ω3 PUFA-enriched diet which lowered this accumulation [START_REF] Javadi | The effect of six different C18 fatty acids on body fat and energy metabolism in mice[END_REF]). Similar results have been obtained in rats, where ω6 PUFA-enriched diet enhanced fat pad weight compared to saturated FA or ω3 PUFA-enriched diets [START_REF] Okuno | Perilla oil prevents the excessive growth of visceral adipose tissue in rats by down-regulating adipocyte differentiation[END_REF]. Interestingly, transgenerational studies developed in rodents have shown that perinatal exposition to high fat diets with high levels of ω6 PUFA induced a progressive accumulation of body fat across generations [START_REF] Massiera | A Western-like fat diet is sufficient to induce a gradual enhancement in fat mass over generations[END_REF]). In the same way, increasing the LA/LNA ratio from 2 to 30 during the preweaning period of guinea pigs led to an increased fat mass at the middle age [START_REF] Pouteau | A low alpha-linolenic intake during early life increases adiposity in the adult guinea pig[END_REF][START_REF] Castaneda-Gutierrez | The guinea pig as a model for metabolic programming of adiposity[END_REF]. In humans, a high LA/LNA PUFA ratio in maternal milk, independently from mother weight, is positively associated with the adiposity of infants at 6 months, 3 and 4 years of age [START_REF] Donahue | Prenatal fatty acid status and child adiposity at age 3 y: results from a US pregnancy cohort[END_REF][START_REF] Moon | Maternal plasma polyunsaturated fatty acid status in late pregnancy is associated with offspring body composition in childhood[END_REF][START_REF] Rudolph | Early infant adipose deposition is positively associated with the n-6 to n-3 fatty acid ratio in human milk independent of maternal BMI[END_REF]. In adults, ω6 PUFA intake correlates positively with the body mass index (BMI) and the associated metabolic syndrome [START_REF] Garaulet | Site-specific differences in the fatty acid composition of abdominal adipose tissue in an obese population from a Mediterranean area: relation with dietary fatty acids, plasma lipid profile, serum insulin, and central obesity[END_REF][START_REF] Williams | Adipose tissue arachidonic acid and the metabolic syndrome in Costa Rican adults[END_REF][START_REF] Inoue | Low serum eicosapentaenoic acid / arachidonic acid ratio in male subjects with visceral obesity[END_REF].

Moreover, in adults both a high intake of ω6 LA and diets associated with a very high LA/LNA PUFA ratios have been associated with the development of several chronic diseases, including cardiovascular, inflammatory and especially overweight/obesity [START_REF] Ailhaud | Fatty acid composition of fats is an early determinant of childhood obesity: a short review and an opinion[END_REF][START_REF] Hibbeln | Healthy intakes of n-3 and n-6 fatty acids: estimations considering worldwide diversity[END_REF], Okuyama, Ichikawa et al. 2007, Okuyama, Ichikawa et al. 2007). Accordingly, new dietary recommendations emphasize higher intake of ω3 PUFAs and lower intake of ω6 PUFAs to improve human health [START_REF] Ailhaud | Temporal changes in dietary fats: role of n-6 polyunsaturated fatty acids in excessive adipose tissue development and relationship to obesity[END_REF][START_REF] Muhlhausler | Omega-6 polyunsaturated fatty acids and the early origins of obesity[END_REF][START_REF] Simopoulos | An Increase in the Omega-6/Omega-3 Fatty Acid Ratio Increases the Risk for Obesity[END_REF][START_REF] Simopoulos | The importance of a balanced omega-6 to omega-3 ratio in the prevention and management of obesity[END_REF] In contrast to ω6 PUFA excessive intake, it has been reported that diets enriched with ω3 PUFAs or ω3 LC-PUFAs lead to beneficial effects and contribute to improve of human health. Cohort analyses have shown normalization of inflammation, lowering plasma triglycerides and amelioration of the cardiovascular metabolic status (Duda, O'Shea et al. 2009, Innes and[START_REF] Innes | The Differential Effects of Eicosapentaenoic Acid and Docosahexaenoic Acid on Cardiometabolic Risk Factors: A Systematic Review[END_REF]. In rodents, inclusion of LNA in LArich high fat diet prevents the enhancement of fat mass [START_REF] Massiera | Arachidonic acid and prostacyclin signaling promote adipose tissue development: a human health concern?[END_REF]. A similar favorable effect is observed for body weight in humans whereas inclusion of EPA and DHA leads to lower cardiovascular disease risk and ameliorates metabolic markers [START_REF] Simopoulos | An Increase in the Omega-6/Omega-3 Fatty Acid Ratio Increases the Risk for Obesity[END_REF]DiNicolantonio 2016, Jayarathne, Koboziev et al. 2017).

Active metabolites of PUFAs

Most biological functions are not attributed directly to PUFAs but rather mediated through some of their oxygenated metabolites. These metabolites, termed oxylipins, are synthesized in a highly regulated manner, mainly by cyclooxygenases (COX-1 and COX-2) and lipoxygenases (3-, 5-, 12-and 15-LOX) activities (Fig. 2), displaying both a high biological potency and a short half-life [START_REF] Funk | Prostaglandins and leukotrienes: advances in eicosanoid biology[END_REF][START_REF] Perez-Chacon | Control of free arachidonic acid levels by phospholipases A2 and lysophospholipid acyltransferases[END_REF]). After synthesis, these oxygenated derivatives can i) undergo new enzymatic modifications, ii) activate soluble intracellular receptors, iii) diffuse out of the cell through the plasma membrane and signal through membrane receptors in a paracrine or autocrine manner, and iv) stored as esterified fatty acids included within cell membranes. In turn, these esterified oxylipins can be hydrolyzed from the membrane under specific stimuli to actively participate into various biological responses [START_REF] Hammond | Esterified eicosanoids: generation, characterization and function[END_REF].

Synthesis of oxylipins takes place from developmental stages to adult life. In mammals, oxylipins can be divided into two main categories (Fig. 2), the eicosanoids (synthesized from C20 PUFAs, i.e. ARA, EPA and dihomo-γ-linolenic acid, DGLA) and the docosanoids (synthesized from C22 PUFAs, i.e. DHA) [START_REF] Massey | Lipidomics of oxidized polyunsaturated fatty acids[END_REF]. A third class of metabolites are derived directly from linoleic acid and composed of two members, 9-and 13-HODEs (hydroxyoctadecadienoic acid). Oxylipins, which differ by their length and double bond configuration, are synthesized by several alternative and consecutive reactions and display crucial signaling functions. Among ARA metabolites, prostacyclin (PGI2) has been shown in mice to trigger in vitro adipose cell differentiation [START_REF] Negrel | Prostacyclin as a potent effector of adipose-cell differentiation[END_REF]) and in vivo adipose tissue formation [START_REF] Massiera | Arachidonic acid and prostacyclin signaling promote adipose tissue development: a human health concern?[END_REF].

There is a direct link between PUFA dietary intake and synthetized oxylipins in plasma and tissues. In rodents, increasing ARA or LA intake increases ω6 PUFAderived oxylipins [START_REF] Pisani | The omega6-fatty acid, arachidonic acid, regulates the conversion of white to brite adipocyte through a prostaglandin/calcium mediated pathway[END_REF][START_REF] Pisani | Disequilibrium of polyunsaturated fatty acids status and its dual effect in modulating adipose tissue development and functions[END_REF][START_REF] Ghandour | Impact of dietary omega3 polyunsaturated fatty acid supplementation on brown and brite adipocyte function[END_REF]), a profile reversed by addition of LNA [START_REF] Ghandour | Impact of dietary omega3 polyunsaturated fatty acid supplementation on brown and brite adipocyte function[END_REF][START_REF] Colson | Diet Supplementation in omega3 Polyunsaturated Fatty Acid Favors an Anti-Inflammatory Basal Environment in Mouse Adipose Tissue[END_REF]. In humans, studies have been mainly focused on ω3-diet supplementation which leads to an increase of ω3-derived metabolites, however the oxylipin levels return to almost normal levels when the ω3-diet supplementation was discontinued [START_REF] Fischer | Dietary omega-3 fatty acids modulate the eicosanoid profile in man primarily via the CYP-epoxygenase pathway[END_REF]. Moreover, in a cohort of hyperlipidemic patients, a 3 months diet supplemented in ω3 PUFAs leads to an increase in ω3-derived metabolites in red blood cell membranes, almost attaining the proportion observed in healthy patients (Harris andVon Schacky 2004, Schmocker, Zhang et al. 2018). This increase in oxylipin synthesis can be acute, as demonstrated recently where both EPA-and DHA-enriched meals induce within 4 hours a potent increase of ω3-derived metabolites in plasma [START_REF] Mcmanus | Differential effects of EPA versus DHA on postprandial vascular function and the plasma oxylipin profile in men[END_REF].

PUFA supplementation and energy expenditure

Few data are available regarding the impact of PUFA and LC-PUFA supplementation on the resting metabolic rate which should be poorly modulated as it is critical for life. In contrast, several studies have analyzed the impact of the different PUFA classes on adaptive thermogenesis and physical activities.

PUFAs and non-shivering thermogenesis

Some studies have been carried out regarding the involvement of PUFA and their metabolites in the activation of BAT and in the conversion of white into brite adipocytes, both in vivo and in vitro [START_REF] Kuda | Omega-3 fatty acids and adipose tissue biology[END_REF]. Recently, our group [START_REF] Pisani | The omega6-fatty acid, arachidonic acid, regulates the conversion of white to brite adipocyte through a prostaglandin/calcium mediated pathway[END_REF][START_REF] Pisani | Disequilibrium of polyunsaturated fatty acids status and its dual effect in modulating adipose tissue development and functions[END_REF][START_REF] Ghandour | Impact of dietary omega3 polyunsaturated fatty acid supplementation on brown and brite adipocyte function[END_REF][START_REF] Colson | Diet Supplementation in omega3 Polyunsaturated Fatty Acid Favors an Anti-Inflammatory Basal Environment in Mouse Adipose Tissue[END_REF]) has shown the importance of the ω6 and ω3 PUFAs and LC-PUFAs dietary intake on adipose tissue development and function. ARA enriched diets, known to favor white adipose tissue formation, prevent the "browning" process to take place in white adipose tissue depots [START_REF] Pisani | The omega6-fatty acid, arachidonic acid, regulates the conversion of white to brite adipocyte through a prostaglandin/calcium mediated pathway[END_REF], suggesting a favorable role of ω3-PUFA supplementation in preventing both processes. Indeed, ω3 PUFA and LC-PUFAs supplementation has beneficial effect on the thermogenic function of adipocytes. A low dietary LA/LNA ratio (4 compared to 25) improved the thermogenic response of BAT and WAT under β3-adrenergic stimulation [START_REF] Ghandour | Impact of dietary omega3 polyunsaturated fatty acid supplementation on brown and brite adipocyte function[END_REF]. The eicosanoid PGF2α appears responsible of these effects, the level of which is influenced by the availability of ω6 and ω3 PUFAs competing at the level of their metabolization.

Few studies showed a direct relationship between ω3 PUFAs enriched diet and induction of brown adipocyte thermogenesis [START_REF] Adkins | Dietary Docosahexaenoic Acid and trans-10, cis-12-Conjugated Linoleic Acid Differentially Alter Oxylipin Profiles in Mouse Periuterine Adipose Tissue[END_REF][START_REF] Kuda | Omega-3 fatty acids and adipose tissue biology[END_REF]). However, it also has been shown that ω3 PUFAs favored the formation of the so-called ''healthy adipocytes" with high number of mitochondria but not expressing UCP1 [START_REF] Kuda | Omega-3 fatty acids and adipose tissue biology[END_REF]. Moreover, other studies have shown the beneficial effect of ω3 PUFA supplementation on glucose metabolism (glucose tolerance, insulin sensitivity) and on BAT activation and/or brite adipocyte formation [START_REF] Kim | Fish oil intake induces UCP1 upregulation in brown and white adipose tissue via the sympathetic nervous system[END_REF][START_REF] Pahlavani | Eicosapentaenoic acid regulates brown adipose tissue metabolism in highfat-fed mice and in clonal brown adipocytes[END_REF]. In vivo studies have shown a global increase in thermogenesis in mice fed fish oil ω3-enriched food, i.e. an increased rectal temperature associated with an increase of brown adipocyte markers in BAT [START_REF] Kim | Eicosapentaenoic Acid Potentiates Brown Thermogenesis through FFAR4-dependent Up-regulation of miR-30b and miR-378[END_REF]. Furthermore, mice fed fish oil (enriched with EPA or DHA) diets also display a strong increase in UCP1 expression in inguinal WAT [START_REF] Kim | Fish oil intake induces UCP1 upregulation in brown and white adipose tissue via the sympathetic nervous system[END_REF]. In this study, was examined the pathway by which ω3 PUFAs could exert their effect using β-adreno-receptor blockers. Indeed, 4 hours of fish oil diet are enough to induce UCP1 in BAT and this effect was abolished when the mice received propranolol as a β-adrenergic receptor blocker [START_REF] Kim | Fish oil intake induces UCP1 upregulation in brown and white adipose tissue via the sympathetic nervous system[END_REF]. In another study, mice fed control diet supplemented with EPA (36 g/kg for 11 weeks) display an increased expression of UCP1 in BAT whereas no effect is observed in WAT [START_REF] Pahlavani | Eicosapentaenoic acid regulates brown adipose tissue metabolism in highfat-fed mice and in clonal brown adipocytes[END_REF]. Although EPA appears important for BAT activity, it cannot be excluded that others ω3 PUFAs present in fish oil are important for this activity.

In addition to their role in brown and brite adipocytes differentiation and function, oxylipins are well-known modulator of white adipocyte development and function. For example, prostacyclin (PGI2) synthesized from ARA can promote adipogenesis. Indeed, secreted prostacyclin from preadipocyte binds to its receptor IP in an autocrine/paracrine way and leads to cAMP pathway activation. This enhances the differentiation of preadipocytes via the expression of C/EBPβ and C/EBPδ (C/Enhancer Binding Proteins) [START_REF] Vassaux | Prostacyclin is a specific effector of adipose cell differentiation. Its dual role as a cAMP-and Ca(2+)-elevating agent[END_REF][START_REF] Massiera | Arachidonic acid and prostacyclin signaling promote adipose tissue development: a human health concern?[END_REF], both of which being critical for the progression of early phase of adipogenesis and the activation of PPARs (Peroxisome Proliferator-Activated Receptors) [START_REF] Brun | Differential activation of adipogenesis by multiple PPAR isoforms[END_REF]. In addition to activating plasma membrane receptors, PGI2 can also signal directly through PPARs activation to promote adipogenesis [START_REF] Forman | 15-Deoxydelta 12, 14-prostaglandin J2 is a ligand for the adipocyte determination factor PPAR gamma[END_REF]. This dual signalling is allowed by the localization of the prostacyclin synthase PGIS at the plasma or at the nuclear membrane [START_REF] Lim | A novel pathway of prostacyclin signaling-hanging out with nuclear receptors[END_REF]. The involvement of PGI2 in adipogenesis have been confirmed in vivo in mice invalidated for its receptor IP. Indeed, these mice displayed less body weight and fat mass compared to wild type mice when fed with a LA-rich diet [START_REF] Massiera | Arachidonic acid and prostacyclin signaling promote adipose tissue development: a human health concern?[END_REF]. In addition to PGI2, some others ARA metabolites generated through cyclooxygenase activities, as PGD2 and 15d-PGJ2, had been described to stimulate adipogenesis [START_REF] Forman | 15-Deoxydelta 12, 14-prostaglandin J2 is a ligand for the adipocyte determination factor PPAR gamma[END_REF].

PUFAs and fever

Infection and inflammation result in a negative energy balance with an increased thermogenesis (also named fever, pyrexia or febrile state) and decreases food intake.

As fever is mediated mainly by leptin, IL-1β, IL-6 and PGE2, all factors modulating production/secretion of these mediators should be able to mediate fever.

IL-1β and IL-6 are two pro-inflammatory cytokines released by immune cells and crucial to the regulation/induction of fever by the brain [START_REF] Kozak | IL-6 and IL-1 beta in fever. Studies using cytokine-deficient (knockout) mice[END_REF]. It is well accepted that ω6 PUFAs leads to the synthesis of inflammatory eicosanoids (including PGE2) and that ω3 PUFAs display anti-inflammatory properties, mainly by the production of anti-inflammatory eicosanoids and docosanoids, but also by the modulation of inflammatory cytokine expression [START_REF] Calder | Polyunsaturated fatty acids and inflammation[END_REF]. Mice feeding diet enriched with fish or rapeseed oil (ω3-enriched oils) displayed a decreased IL-1β and/or IL-6 levels in blood [START_REF] Sadeghi | Dietary lipids modify the cytokine response to bacterial lipopolysaccharide in mice[END_REF] as well as in brain [START_REF] Delpech | Dietary n-3 PUFAs Deficiency Increases Vulnerability to Inflammation-Induced Spatial Memory Impairment[END_REF] under endotoxemia conditions. This effect of ω3 PUFAs still under debate as it was not already found under healthy conditions [START_REF] Colson | Diet Supplementation in omega3 Polyunsaturated Fatty Acid Favors an Anti-Inflammatory Basal Environment in Mouse Adipose Tissue[END_REF]) but mainly under pathological or inflammatory conditions.

PGE2 is an indispensable brain mediator of fever as demonstrated in studies with targeted disruption of genes encoding either PGE2-synthesizing enzymes [START_REF] Engblom | Microsomal prostaglandin E synthase-1 is the central switch during immune-induced pyresis[END_REF] or PGE2 receptors [START_REF] Ushikubi | Impaired febrile response in mice lacking the prostaglandin E receptor subtype EP3[END_REF], or using pharmaceutic compounds such as aspirin, salicylate, indomethacin or celecoxib [START_REF] Romanovsky | do we have a new target for anti-inflammatory and antipyretic therapy?[END_REF]. In addition, others prostaglandin and non-prostaglandin metabolites derived from arachidonic acid have been involved in the central control of fever in response to inflammation [START_REF] Kozak | Non-prostaglandin eicosanoids in fever and anapyrexia[END_REF]. Due to the tight regulation of oxylipin synthesis in face of the quantity and quality of dietary PUFAs, a direct correlation exists between ω6/ω3 PUFAs and LC-PUFAs ratio and the various pathways involving oxylipins in fever.

Leptin has been characterized in murine models as a fever mediator of body temperature [START_REF] Luheshi | Leptin actions on food intake and body temperature are mediated by IL-1[END_REF]. Little is known about the detailed mechanism linking leptin and pyrexia. Leptin acts at the level of hypothalamus, maybe by modulation of IL-1β expression and IL-6 [START_REF] Pohl | Leptin modulates the late fever response to LPS in diet-induced obese animals[END_REF], which induces a decrease in heat loss and thus an increase in body temperature [START_REF] Fischer | Leptin Raises Defended Body Temperature without Activating Thermogenesis[END_REF][START_REF] Kaiyala | Physiological role for leptin in the control of thermal conductance[END_REF]. As PUFA intake modulated fat mass development which is directly correlated to leptin secretion [START_REF] Garaulet | Site-specific differences in the fatty acid composition of abdominal adipose tissue in an obese population from a Mediterranean area: relation with dietary fatty acids, plasma lipid profile, serum insulin, and central obesity[END_REF][START_REF] Williams | Adipose tissue arachidonic acid and the metabolic syndrome in Costa Rican adults[END_REF][START_REF] Inoue | Low serum eicosapentaenoic acid / arachidonic acid ratio in male subjects with visceral obesity[END_REF], we assume that regulation of fat mass by dietary PUFAs indirectly modulated fever and thus energy expense during inflammation/infection. Interestingly, DHA has been characterized as a determinant of leptin gene expression [START_REF] Shen | Epigenetic modification of the leptin promoter in diet-induced obese mice and the effects of N-3 polyunsaturated fatty acids[END_REF], but human studies failed to clearly correlate ω3 PUFA intake to leptin circulating level in diabetic and healthy populations [START_REF] Kratz | The impact of dietary fat composition on serum leptin concentrations in healthy nonobese men and women[END_REF][START_REF] Stirban | Effects of n-3 polyunsaturated fatty acids (PUFAs) on circulating adiponectin and leptin in subjects with type 2 diabetes mellitus[END_REF].

PUFAs and muscle function

Physical activity increases with body weight as body mass to move is increased. Indeed, even if amplitude and velocity are lowered in obese and overweight subjects, the whole strength developed is similar to lean people [START_REF] Lazzer | Assessment of energy expenditure associated with physical activities in free-living obese and nonobese adolescents[END_REF].

Obesity is still associated (as a cause or a consequence) with a sedentary lifestyle. This behavior, involving a chronic decrease of physical activity frequency, drives de facto a decrease in energy expenditure associated to limited daily life activities [START_REF] Carneiro | Is Obesity Associated with Altered Energy Expenditure?[END_REF]. Therefore, as PUFAs intake has been related to fat mass development, if causality between decreased activities and PUFAs intake exists; it should be rather indirect.

Nevertheless, several studies have described a relationship between dietary PUFAs and muscle mass and function, independent of body weight change [START_REF] Jeromson | Omega-3 Fatty Acids and Skeletal Muscle Health[END_REF][START_REF] Gammone | Omega-3 Polyunsaturated Fatty Acids: Benefits and Endpoints in Sport[END_REF]. Firstly, PUFAs seem to affect muscle anabolism, including muscle growth and regeneration [START_REF] Tachtsis | Potential Roles of n-3 PUFAs during Skeletal Muscle Growth and Regeneration[END_REF]. Other studies have shown an increased muscle protein synthesis in both young and older subjects undergoing a daily administration of ω3 PUFAs for eight weeks [START_REF] Smith | Omega-3 polyunsaturated fatty acids augment the muscle protein anabolic response to hyperinsulinaemia-hyperaminoacidaemia in healthy young and middle-aged men and women[END_REF]. A similar result on muscle mass is obtained after six months of ω3 PUFAs supplementation in older people with an additional increase in muscle strength [START_REF] Smith | Fish oilderived n-3 PUFA therapy increases muscle mass and function in healthy older adults[END_REF]. Secondly, ω3 PUFA intake ameliorates muscle metabolism and function. In mice submitted to a high fat diet, a four-week daily supplementation with ω3 PUFAs reversed lipid accretion and increased oxidative metabolism in muscle. Interestingly, in this study the authors describe an increase in BAT mass which can be relied to the positive effect of ω3 PUFA intake [START_REF] Philp | Dietary enrichment with fish oil prevents high fat-induced metabolic dysfunction in skeletal muscle in mice[END_REF].

Overall, even if ω6 and ω3 PUFA intakes do not modulate physical activities per se, they have an impact on skeletal muscle capacity and can influence indirectly energy expenditure during exercise.

Conclusions

ω6 and ω3 PUFAs and LC-PUFAs can be considered as two important actors regulating energy expenditure. Neither exclusive positive nor unfavorable involvement of the different PUFAs is observed, as each of them can promote or inhibit different kinds of energy expenditure (Fig. 3). Nevertheless, excessive dietary intake of ω6 PUFAs can overall be considered as negative regulators and adequate intake of ω3 PUFAs as positive counter-regulators of various metabolic events. Indeed, excess of ω6 PUFAs or insufficient ω3 PUFAs intake favor fat mass development which can be relied to a decrease in both non-shivering thermogenesis and physical activity due to sedentary lifestyle. Conversely, ω3 PUFAs supplementation is limiting fat mass increase as well as promoting thermogenesis and physical activity. Of note, as both ω6 and ω3 PUFAs are essential and metabolically related to each other (Fig. 2), it is difficult to bring unequivocal and direct effect of ω6 PUFAs distinct from ω3 PUFAs on the different types of energy expenditure (Fig. 3). Another critical issue of the reported studies has been the use of high ω6 PUFAs intake and/or various ω3 PUFAs supplementation which appear widely different from nutritional recommendations.

  ) acids. These PUFAs (LA and LNA) and long-chain PUFAs (LC-PUFAs) (ARA, DGLA, EPA and DHA) trigger a variety of biological responses and are required for normal development.

  Overall, in 2010, the French Agency for Food, Environmental and Occupational Health & Safety (ANSES) recommended an intake of 4% of LA, 1% of LNA and 250 //www.anses.fr/en/content/fats), this nutritional recommendation appearing adequate to regulate an optimal impact on energy storage and expenditure.

Figure 1 .

 1 Figure legends.

Figure 2 .

 2 Figure 2. Major pathways of PUFA metabolization, from dietary fatty acids to major members of the oxylipin family, including the commons enzymatic hubs of metabolization between ω6 and ω3 PUFAs. COX: cyclooxygenase; LOX: lipoxygenase; HODEs: hydroxyoctadecadienoic acids; LA: linoleic acid; LNA: alphalinolenic acid; ARA: arachidonic acid; EPA: eicosapentaenoic acid; DHA: docosahexaenoic acid. The number of carbon atoms is indicated.

Figure 3 .

 3 Figure 3. Major thermogenic tissues, interrelationships and modulation by ω6 and ω3 PUFAs and LC-PUFAs. Red rounded-head arrow: negative action. Green arrow: positive action. Black arrow: modulatory pathways between organs. Grey arrow: environmental effectors. Flame: heat production.


Claude Pasquier
email: claude.pasquier@univ-cotedazur.fr

Vincent Guerlais

Denis Pallez

Raphaël Rapetti-Mauss

Olivier Soriani

Identification of active modules in interaction networks using node2vec network embedding

The identification of condition-specific gene sets from transcriptomic experiments is important to reveal regulatory and signaling mechanisms associated with a given cellular response. Statistical approaches using only expression data allow the identification of genes whose expression is most altered between different conditions. However, a phenotype is rarely a direct consequence of the activity of a single gene, but rather reflects the interplay of several genes to carry out certain molecular processes. Many methods have been proposed to analyze the activity of genes in light of our knowledge of their molecular interactions. However, existing methods have many limitations that make them of limited use to biologists: they detect modules that are too large, too small, or they require the users to specify a priori the size of the modules they are looking for.

We propose AMINE (Active Module Identification through Network Embedding), an efficient method for the identification of active modules. Experiments carried out on artificial data sets show that the results obtained are more reliable than many available methods. Moreover, the size of the modules to be identified is not a fixed parameter of the method and does not need to be specified; rather, it adjusts according to the size of the modules to be found. The applications carried out on real datasets show that the method enables to find important genes already highlighted by approaches solely based on gene variations, but also to identify new groups of genes of high interest. In addition, AMINE method can be used as a web service on your own data (http://amine.i3s.unice.fr).

Introduction

Current high-throughput technologies are now capable of reliably quantifying, at the scale of an entire organism, the molecular changes that arise in response to diseases or environmental disturbances. In order to extract from this pool of data the genes most related to the process under study (and therefore of most interest to the biologist), statistical methods are generally used to associate genes with numerical values reflecting the extent of their variation. In most studies, the genes considered of most interest are the ones whose relative differences in expression, or fold changes, are the largest. Unfortunately, the raw fold change is unreliable because it does not take into account the uncertainty inherent in gene expression measurements. To overcome this uncertainty, existing methods calculate a p-value to reflect the statistical significance of the variation.

Selecting the genes of interest on the basis of fold changes, p-values or a combination of both makes it possible to compile a list of genes whose expression varies most significantly. However, this procedure fails to identify genes whose combined action is essential in the process under study but whose individual scores are too low.

Though, as pinpointed by [START_REF] Rapaport | Classification of microarray data using gene networks[END_REF], "a small but coherent difference in the expression of all the genes in a pathway should be more significant than a larger difference occurring in unrelated genes". Arising from this observation, many methods have been proposed to analyze gene activity in the light of our knowledge about their molecular interactions. These pertinent subnetwork are named "context-dependent active subnetworks" [START_REF] He | Comparison of statistical methods for subnetwork detection in the integration of gene expression and protein interaction network[END_REF], "functional module" [START_REF] Beisser | BioNet: an R-Package for the functional analysis of biological networks[END_REF], "maximal scoring subgraph" [START_REF] Dittrich | Identifying functional modules in protein-protein interaction networks: an integrated exact approach[END_REF] or "altered subnetworks" [START_REF] Reyna | Hierarchical HotNet: identifying hierarchies of altered subnetworks[END_REF]. The underlying idea is to identify a pertinent module of genes by simultaneously taking into account two criteria: one based on a measurement of genes activity and the other one reflecting the proximity between the genes in the module. One of the challenges is to define an appropriate scoring strategy based on these two criteria. [START_REF] Nguyen | A comprehensive survey of tools and software for active subnetwork identification[END_REF] classified main computational methods for solving the active subnetworks identification problem in six categories: (i) greedy algorithms, (ii) random walk algorithms, (iii) diffusion emulation models, (iv) evolutionary algorithms, (v) maximal clique1 identification and (vi) clustering based methods. First two methods are simple and rapid but are highly dependent to the starting point of the algorithm that does not guarantee to reach global optima. Conversely, methods (iii) and (iv) are able to find global optima (in accordance with the scoring system used) or an approximation of it at the prize of a computational burden. Methods (v) do not fully answer the initial issue as it is probably not true that each gene involved in a biological process interacts with all the others. Finally, methods (vi) offers the advantage of being based on existing clustering algorithms, but these rely on a distance between objects that must be determined. Moreover, most of the clustering algorithms require to determine a priori the number of clusters to build, which is challenging. The commonality between all these methods is that their effectiveness is very dependent of the network topology. Unfortunately, it is known that molecular interaction networks are noisy and incomplete [START_REF] De Las Rivas | Protein-protein interactions essentials: key concepts to building and analyzing interactome networks[END_REF]. In recent years, network embedding [START_REF] Cui | A survey on network embedding[END_REF] has proven to be a powerful network analysis approach by generating a very informative and compact vector representation for each vertex v in the network. The approach was initially considered as part of dimensionality reduction techniques (reducing for example a |v|×|v| adjacency matrix into a |v|× m matrix where m≪|v|). This dimensionality reduction allows to reduce noise and map nodes in a vector space in which distances between nodes accurately reflect their proximity in the original network. One of the most representative technique for network embedding is Node2vec [START_REF] Grover | node2vec: Scalable feature learning for networks[END_REF].

Recent advances on deep learning has led to a plethora of methods based on deep neural networks for learning graph representations, methods that are often inspired by the learning of word embedding [START_REF] Mikolov | Distributed representations of words and phrases and their compositionality[END_REF]. Works on word embedding can be seen as learning linear sequences (word sequences). It has been shown that the resulting compact vector representations are capable of capturing rich semantic information about natural language. Processing graph structures is much more complicated. A popular approach is to convert a complex graph structure with a rich topology into a set of linear structures and then use a word embedding method to calculate the vector representation of each node.

As mentioned above, the identification of active modules requires the simultaneous consideration of two criteria. In existing methods, measurements of gene activity and their network proximity are either combined to form a single metric or optimized simultaneously using multiobjective algorithms [START_REF] Corrêa | Population-based meta-heuristic for active modules identification[END_REF]. When working on embedded networks, the proximity facet is embedded in the vector space. It is then possible to focus on the detection of subspaces containing genes that have a high activity. Following this line, we propose AMINE (Active Module Identification through Network Embedding), a new and efficient method for active module detection based on Node2vec [START_REF] Grover | node2vec: Scalable feature learning for networks[END_REF]. Our method uses a greedy approach to build increasingly large clusters of nodes based on the similarity of their encoding vectors and to evaluate them according to a metric taking into account the activity of the contained nodes.

On artificially generated datasets, the effectiveness of the method compares well with several recent algorithms. On real datasets, AMINE allows to complement the results obtained with classical approaches by identifying new groups of genes of great interest.

Results

Evaluation of AMINE on artificial data generated by Robinson et al. (2017)

Many studies dealing with the identification of active modules have tested their methods on datasets generated by themselves and which are, at times, difficult to reproduce. [START_REF] Robinson | Incorporating interaction networks into the determination of functionally related hit genes in genomic experiments with Markov random fields[END_REF] are among the few to give access to all materials used to test the MRF method they proposed. These materials contain the graph itself, the p-values associated with the nodes and the modules to be identified. It gives us the opportunity to apply our method on exactly the same data.

The simulated experiment used to evaluate the MRF method [START_REF] Robinson | Incorporating interaction networks into the determination of functionally related hit genes in genomic experiments with Markov random fields[END_REF] consists of a set of 1000 scale free graphs, each containing 1000 vertices associated with values simulated from a standard uniform distribution. In this dataset, there are three modules to discover (the hit modules), each containing 10 vertices. To simulate the fact that the vertices belonging to these modules represent differentially expressed genes, and are therefore associated with low p-values, these vertices are assigned simulated values from a truncated Gaussian distribution with mean 0 and standard deviation equals to 10 -6 . [START_REF] Robinson | Incorporating interaction networks into the determination of functionally related hit genes in genomic experiments with Markov random fields[END_REF] compared their MRF method to NePhe [START_REF] Wang | A network-based integrative approach to prioritize reliable hits from multiple genome-wide RNAi screens in Drosophila[END_REF], Knode [START_REF] Cornish | SANTA: quantifying the functional content of molecular networks[END_REF] and BioNet [START_REF] Beisser | BioNet: an R-Package for the functional analysis of biological networks[END_REF] and report that MRF gives the best performance in term of recall. Since it is known that there are exactly 30 true hits in the dataset, the authors rate the different methods by considering only the proportion of true hits in each hit list of size 30.

We ran AMINE on these data to detect the three most significant modules. Since our method does not consider the size of expected modules as an input parameter, it can identify modules of various sizes. When considering only the recall criterion (proportion of hit nodes identified) as Robinson et al. did, AMINE is disadvantaged when it identifies small modules and favored when it identifies large modules. In order to have a fair comparison between AMINE and MRF, we apply the following rule: if the number of genes belonging to the three identified modules is greater than 30, then only 30 genes randomly chosen are considered. In the case where the three identified modules contain less than 30 genes, the recall obtained by our method is underestimated. Despite this disadvantage, the predictions made by AMINE obtained a slightly better recall than the one obtained by MRF (Fig. 1A, last boxplot ; the detail of the evaluation performed by Robinson et al. can be found in Fig. 7 of their article [START_REF] Robinson | Incorporating interaction networks into the determination of functionally related hit genes in genomic experiments with Markov random fields[END_REF]). The recall without truncation (number of true positive divided by the number of hit nodes), the precision (number of true positive divided by the number of nodes identified) and the F1 score (harmonic mean of precision and recall) were also plotted on the same figure. The median value of the F1 score is 0.76 while the minimum and maximum values are 0.31 and 0.94 respectively. The total number of genes in the three identified modules range from 18 to 46 with a median of 28 (Fig. 1B). This means that AMINE is able to identify modules close to the ground truth (although slightly smaller) without the need to specify their size a priori.

Validation of the method on artificial dense networks

It has already been shown that simulating a biological network is a very difficult task [START_REF] Pavlopoulos | Using graph theory to analyze biological networks[END_REF]. However, we found that generating artificial graphs using Barabasi-Albert model of preferential attachment [START_REF] Barabási | Emergence of scaling in random networks[END_REF], as performed in different studies (for example, the articles of [START_REF] Cornish | SANTA: quantifying the functional content of molecular networks[END_REF] and [START_REF] Robinson | Incorporating interaction networks into the determination of functionally related hit genes in genomic experiments with Markov random fields[END_REF]) is too far from a real interaction graph for the results to be extrapolated (see supplementary Fig. S1 for an example of such sparse graph).

Using an extended version of the Barabási-Albert model of preferential attachment [START_REF] Albert | Topology of evolving networks: local events and universality[END_REF], we generated several artificial networks with topologies relatively close to real interaction networks and one hit module to discover with size of 10 or 20 nodes (see materials and methods for details relative to the generation of dense networks and supplementary Fig. S2 for an example of a generated dense network). The performance of AMINE was compared with the methods COSINE [START_REF] Ma | COSINE: COndition-SpecIfic sub-NEtwork identification using a global optimization method[END_REF], BioNet [START_REF] Beisser | BioNet: an R-Package for the functional analysis of biological networks[END_REF], GiGA [START_REF] Breitling | Graph-based iterative Group Analysis enhances microarray interpretation[END_REF]) and a baseline consisting in simply picking the genes with the lowest p-values. The results are shown in Fig. 2 for networks comprising 1000 vertices and modules of size 10 and in Fig. 3 for networks with 1000 vertices and modules of size 20.

We can observe that identifying module on denser network is a process much more complicated as the median F1 score drops from 0.76 on a sparse graph to values just above 0.5, depending of the size of the module. The scores of the other methods are much lower in the case of modules of size 10 ; some even score lower than the baseline. In the case of COSINE and BioNet, these poor results come with the prediction of large modules whose size exceeds 100 genes (Fig. 2B and Fig. 3B). For modules of size 20, the F1 scores of GiGA and AMINE are very close (Fig. 3B). However, it should be noted that GiGA uses a parameter that determines the maximum size of the module to be identified. In our experiment, we have set the maximum size to the size of the module to be identified, which, of course, helps the procedure. As shown in Fig. 2C and Fig. 3C, for GiGA, the size of the identified modules is always below the maximum. For AMINE, without any indication on the size of the modules being searched, we can see that the method predicts modules with a median size close to the ground truth size. This is a significant advantage of our method since biologists do not know in advance the size of the modules they are looking for.

Scalability of the method

To test to what extent our method is able to scale to larger networks, we applied it to an artificial network of 10,000 vertices, that was generated using the same parameters previously defined. The processing time increases from 1 minute for a dense network of 1,000 nodes to 30 minutes for a dense network of 10,000 nodes. This processing time is acceptable as real biological networks are close to this size.

The distribution of the F1 scores obtained by AMINE as well as the distribution of the size of the identified modules are shown in supplementary Fig. S3-S4. For modules of size 10, the median F1 is in agreement with that obtained on networks of size 1,000. The overall performance is however less good because it is penalized by the fact that, in many cases, the module is completely missed (which explains why the second quartile starts at zero in the leftmost boxplot of supplementary Fig. S3). This can be explained by the fact that, as the number of nodes in the network increases, the number of non hit nodes that are associated with a random value higher than the values associated with hit genes increases. So, the probability that random modules score higher than the hit module increases too. This effect is less important for modules of size 20, for which the method works as well as for networks of 1,000 nodes. Regarding sizes, they are still close to the ground truth, although there is a greater spread of values for modules of size 20.

Validation using a real gene expression dataset

In order to test the ability of AMINE to identify relevant biological functions, we downloaded from Gene Expression Omnibus a dataset relative to a study aimed at characterizing processes and genes associated to metastatic spreading in pancreatic ductal adenocarcinoma (PDAC). With a 5-year survival that has not significantly evolved for 30 years despite progresses in anticancer therapies (<6%), PDAC is a cancer with one of the bleakest prognoses of the most fatal cancers. In the study of [START_REF] Chiou | BLIMP1 induces transient metastatic heterogeneity in pancreatic cancer[END_REF], the authors compared two populations of primary PDAC cells according to the expression of HMGA2, a gene associated to a high metastatic potency and poor outcome in several cancers, including PDAC. Gene expression AMINE profiling was carried out from six pairs of HMGA2+/HMGA2-cell populations, each pair originating from PDAC primary tumors spontaneously generated in a genetically engineered PDAC mice model (PKC mice). The profiling generated 193 statically relevant modules (with associated p-values < 0.05) and we decided to focus on the five ones figuring at the top of the list to verify whether they were relevant to metastatic process in PDAC from an ontological point of view. The functions and signaling pathways associated to these five modules were explored using STRING database [START_REF] Szklarczyk | STRING v11: protein-protein association networks with increased coverage, supporting functional discovery in genome-wide experimental datasets[END_REF] and are shown in Table 1 and Fig. 4.

Extra cellular Matrix (ECM) organization and ECM cell interaction

A hallmark of PDAC is a pronounced collagen-rich fibrotic extracellular matrix (ECM) produced by fibroblasts and cancer cells, known as the desmoplastic reaction. The neoplastic epithelium exists within a dense stroma, which is recognized as a critical mediator of disease progression through direct effects of ECM on cancer cells [START_REF] Hosein | Pancreatic cancer stroma: an update on therapeutic targeting strategies[END_REF]. Interestingly, three out of the five modules generated by AMINE were associated to ECM matrix organization (module 1 and 4) and cell interaction with ECM (module 3). Module 1 was more specifically linked to collagen fibril organization, one of the major constituent of PDAC ECM. Indeed, Collagen contributes to tumor cell aggressiveness, metastatic process and chemoresistance [START_REF] Shields | Biochemical role of the collagen-rich tumour microenvironment in pancreatic cancer progression[END_REF][START_REF] Hessmann | Microenvironmental determinants of pancreatic cancer[END_REF]. Interestingly, Module 3 brings together genes involved in the regulation of cancer cell interaction with ECM through focal adhesion kinases and PI3K/AKT pathways (Fig. 4 and Table 1). Based on the available literature, these processes have been strongly involved in the aggressiveness of PDAC cell and the development of metastasis [START_REF] Jiang | Targeting focal adhesion kinase renders pancreatic cancers responsive to checkpoint immunotherapy[END_REF][START_REF] Jiang | Development of resistance to FAK inhibition in pancreatic cancer is linked to stromal depletion[END_REF].

Response to Hypoxia

Extensive desmoplasia and hypo vascularization within PDAC results in significant intra-tumoral hypoxia (low oxygen) that contributes to its aggressiveness, therapeutic resistance, and high mortality [START_REF] Hollinshead | Respiratory supercomplexes promote mitochondrial efficiency and growth in severely hypoxic pancreatic cancer[END_REF][START_REF] Koong | Pancreatic tumors show high levels of hypoxia[END_REF]. Functional enrichment of modules 2 and 5 raised hypoxia-triggered functions, i.e. VEGF-and HIF-1-dependent pathways. These pathways drive angiogenesis, metabolism adaptation of cancer cell to hypoxia (Warburgh effect), cell cycle inhibition, enhanced migration and metastatic progression (Table 1 and Fig. 4). These results are in good agreement with the literature on PDAC; for example, these pathways are over-represented in genome-wide transcriptome profiling from ex-vivo human PDAC [START_REF] Ghaderi | Genome-wide transcriptome profiling of ex-vivo precision-cut slices from human pancreatic ductal adenocarcinoma[END_REF][START_REF] Shah | Hypoxia: Friend or Foe for drug delivery in Pancreatic Cancer[END_REF]. More interestingly, hypoxia-, VEGF-and HIF-1-associated pathways were repressed in PDAC cells in which BLIMP1, one of the most overexpressed gene in HMGA2+ cell subpopulation, was silenced [START_REF] Chiou | BLIMP1 induces transient metastatic heterogeneity in pancreatic cancer[END_REF].

Altogether, these results validate our methods since non-oriented analysis of genes deregulated in pro-metastatic PDAC cells by AMINE retrieves gene modules involved in highly relevant functions in the context of the disease.

In vitro functional validation of Blimp1-associated module in human PDAC cells

Blimp1 is one of the most overexpressed gene in pro-metastatic HMGA2+ PDAC cells [START_REF] Chiou | BLIMP1 induces transient metastatic heterogeneity in pancreatic cancer[END_REF]. In their study, Chiou et al. analyzed the consequences of Blimp1 silencing in mice PDAC cells to unveil its function in disease progression. Based on in silico, in vitro and in vivo experiments, the authors concluded that Blimp1 acts as a driver of the metastatic ability of PDAC cells. In particular, they found that Blimp1 is a Hypoxia/Hif-regulated gene in human and murine PDAC which is in a good agreement with functions recovered in modules 2 and 5 raised by AMINE processing (Table 1 and Fig. 4). Surprisingly, Blimp1 was not included in these modules, but was indeed detected in a module of 4 genes (module number 169 with a p-value of 0.044; Fig. 5A). Functional enrichment of this module unveiled functions associated to immune response, including regulation inflammation, interleukin production and Th17 cell differentiation (Table 2 and Fig. 5A). Interestingly, it is known that neoantigen expression in pancreatic ductal adenocarcinoma (PDAC) results in exacerbation of an inflammatory microenvironment that drives disease progression and metastasis [START_REF] Hegde | Dendritic cell paucity leads to dysfunctional immune surveillance in pancreatic cancer[END_REF]. It was therefore tempting to validate this result using a series of functional experiments. In this perspective, we first explored a RNA-Seq experiment performed by Chiou et al. revealing deregulated genes in Blimp1-silenced PDAC cells compared to control [START_REF] Chiou | BLIMP1 induces transient metastatic heterogeneity in pancreatic cancer[END_REF]. AMINE profiling of genes negatively regulated by Blimp1 silencing revealed 345 modules with associated p-values < 0.05. Among the ten best modules, we found that modules 2 (p-value < 1.01e-11) and 9 (p-value < 2.83e-8) were associated to cytokine production and inflammatory process (Table 3 and Fig. 5B) after functional enrichment. Next, in order to confirm the putative involvement of Blimp1 in epithelial cancer cell inflammatory process in vitro, we silenced BLIMP1 in MIA PaCa-2 cells, a human PDAC cell line, using siRNA silencing, (Fig. 6A), and explored how it modified the profile of cytokine secretion using a cytokine profiling array. Indeed, we found that Blimp1 repression triggered the production of IL-18Bpa and angiogenin, two anti-inflammatory factors [START_REF] Lee | Angiogenin reduces immune inflammation via inhibition of TANK-binding kinase 1 expression in human corneal fibroblast cells[END_REF] and reduced the secretion IL-6, a major proinflammatory interleukin [START_REF] Tanaka | IL-6 in inflammation, immunity, and disease[END_REF] (Fig. 6B).

Altogether, these results further validate AMINE as valuable method to detect relevant functional modules from large experimental datasets, and unveil a new function of Blimp1 in PDAC-related inflammatory process.

Discussion

This paper proposes a new method for identifying gene modules that are activated as a result of a state shift caused by a biological experiment. Our method, called AMINE, uses as inputs, on the one hand, the ultimate result of any RNA-Seq analysis pipeline which is the differential expression of genes, and on the other hand, a network modeling the interactions between genes. Regarding the second point, the version of our method that is executable (at the address http://amine.i3s.unice.fr/) makes use of the protein interaction data provided by the STRING database for four model organisms: Caenorhabditis elegans, Drosophila melanogaster, Homo sapiens and Mus musculus. Thus, a user only needs to provide the differential gene expression data generated by the pipeline of his choice. From a very simple interface (supplementary Fig. S5), he only has to specify the name of the organism analyzed, the file on which the data are located and the Id of the columns containing the genes' names, the p-values and optionally the fold changes to be able to launch the process. The address of the page containing the results is e-mailed to the user when the processing is completed. On the result page, the most significant modules are listed, however, all the modules found can be downloaded as an Excel document consisting of two sheets. The first sheet, named "list of modules" contains the list of all modules found. The results are presented in 4 columns containing the module number, the list of genes in the module, the s score of the module and the associated p-value. The second sheet, named "genes to modules", is composed of two columns: the first one contains the name of a gene and the second one, the module to which it belongs.

Despite the large number of methods developed over the last 20 years, AMINE identifies, with a higher accuracy than its competitors, modules created computationally on datasets intended to mimic the topology of biological networks. Extrapolating these results to a measure of accuracy on real datasets is very difficult. There is no method to ensure that good predictions on artificial data translate into good predictions on real datasets. However, we have made a special effort to ensure that our simulations are close to real datasets. The networks we generate, with the parameters presented in this article, are closer to a real interaction network than the networks used by some competing methods.

However, these efforts run up against reality. Indeed, it is known that the interaction networks stored in public databases are both incomplete and contain erroneous interactions. In their study, Von [START_REF] Mering | Comparative assessment of large-scale data sets of protein-protein interactions[END_REF] estimate that, for Saccharomyces cerevisiae, the protein-protein interaction data (PPI) reported in public databases account for only one third of existing interactions. This observation leads us to believe that methods based very precisely on the topology of networks are not to their advantage. We think first of all of the methods based on the identification of cliques. Their performance is more than questionable if we consider that a large part of the interactions between proteins are unknown. If we consider that the missing interactions are randomly distributed on the graph, we can estimate that all the paths on the graph are impacted in the same way and thus that the methods based on random walks could be the least affected. Intuitively, we can indeed argue that, in a graph on which a certain proportion of the edges have been randomly deleted, if, from a source node A, random walks allow on average to reach node B before node C, then, on the complete graph, node B will probably always be closer to node A than node C. The other problem with methods based on graph traversals is that PPI networks are "smallworld" networks, meaning that the neighbors of a given node are likely to be neighbors of each other, and most nodes can be reached from every other node by a small number of hops. Thus, any method that relies on graph traversals will find that a large portion of the network is close to any typical node [START_REF] Cao | Going the distance for protein function prediction: a new distance metric for protein interaction networks[END_REF].

Based on these considerations, we argue that network embedding methods can provide the backbone of a more reliable method by estimating distances between nodes that take into account the entire topology of the graph and, moreover, are little affected by the proportion of missing edges. Our method works on an embedding of an interaction network by adopting a greedy algorithm and an active subnetwork relevance measure defined in other papers [START_REF] Ideker | Discovering regulatory and signalling circuits in molecular interaction networks[END_REF]. The great advantage of our method is that it does not require any parameterization; it is not even necessary to indicate the number of modules to be identified or the size of the modules.

We have checked that our method performs well on artificial datasets and compares favorably with existing methods that are the current state of the art. We then processed a real dataset from a study focused on PDAC, on which AMINE retrieved modules associated with functions involved in PDAC metastatic process. Moreover, we highlighted functions that could not be detected using traditional approaches consisting in analyzing only the most differentially expressed genes. Our studies show that AMINE can identify modules corresponding to functions not revealed by traditional approaches. Indeed, we found that Blimp1, one of the most up-regulated gene in highly metastatic cells, participated in the regulation of pro-inflammatory process, a result confirmed by in vitro by the silencing of Blimp1 in human PDAC cells. However, we stress that our method is not an alternative to methods based on the identification of the most differentially expressed genes, but rather a complement to these approaches.

Materials and Methods

The method AMINE predicts active modules from data consisting of a background knowledge about gene interactions and measurements representing, in the specific context of a given experiment, indicators of the involvement of genes in the studied process. This concept of geneinvolvement is materialized by a p-value which quantifies, for each gene, the statistical significance of its variation (Fig. 7A).

Data about gene interactions and gene variations are merged to generate an attributed gene network in which genes are annotated with a numeric attribute representing the extent of their variation (Fig. 7B). Mathematically, the dataset is represented as an attributed graph G=(V , E , λ) consisting of a set of vertices V (also called node, that symbolize the genes), a set of edges E ⊆ {(u , v ) ∈V 2 ∨u ≠ v } and a value function λ ( v ):V → R which associates a value p ∈ R to each vertex v ∈ V . An induced subgraph of G is a subset S of the vertices of G together with those edges of G with both endpoints in S. Many active module detection algorithms focus on identifying induced subgraphs whose values associated with their nodes stand out from the values associated with the other nodes of the graph. We hypothesize that focusing heavily on the detection of connected sets of genes may not be optimal, given the fact that the interactions between genes described in the databases are still largely incomplete. For this, we adopt a definition of a module that is closer to the one used in cluster analysis: object that are grouped together (in a module) are more similar to each other than to those in other groups. The notion of similarity encompasses a component taking into account the distance on the graph between the vertices belonging to a module and a significant nearness between the values associated to these vertices.

Scoring of a subgraph

Let p i = λ ( v i ) be the associated p-value of vertex v i . We aggregate the p-values associated to the nodes of a subgraph with Stouffer's Z method; the same strategy used by [START_REF] Ideker | Discovering regulatory and signalling circuits in molecular interaction networks[END_REF]. If we

let z ( v i )=Φ -1 ( 1 -λ ( v i ))
, where Φ is the standard normal cumulative distribution function, then, the aggregate z-score z a (G' ) for an induced subgraph G ' ⊆G composed of k vertices, is computed with:

z a (G' )= 1 √k ∑ v i ∈G ' z ( v i )
In order to get a subnetwork which has higher aggregation z-score compared with a random set of vertices, we define, still following the same methodology as [START_REF] Ideker | Discovering regulatory and signalling circuits in molecular interaction networks[END_REF], a corrected score s (G' ) of a subgraph G ' with:

s (G' )= z a (G ') -μ k σ k
where the mean μ k and standard deviation σ k are computed based on a Monte-Carlo approach, taking 10,000 rounds of randomly sampling a connected subgraph of k vertices from V . From s (G' ) , we can easily compute the probability of observing in G, a subnetwork of the same size as G ' with a corrected score at least as extreme as the one observed. This is given by the the one-sided p-value:

pvalue(G ' )=1-Φ ( s (G' )).

Network embedding

A network embedding method is a function ψ :V → R m that associate to each vertex v of the graph a vector d of size m. Node2vec (Grover and Leskovec, 2016) uses a biased random walk procedure which efficiently generates diverse neighborhoods of a given node. Node contexts are then processed with the word2vec method [START_REF] Mikolov | Efficient estimation of word representations in vector space[END_REF]. Node2Vec uses two parameters to control the walks. Intuitively, these parameters control how the walk explores and leaves the neighborhood of starting nodes. They allow a tuning between outward exploration and local walking.

However, in our case study, it may indeed be interesting, instead of using biases that only considered the topology of the network, to use the data associated to nodes, i.e. the value of p. The idea is to bias the walk so that, when the walker is located on a node, transitions to nodes with similar values of p (Fig. 7C) are favored. As p represents a p-value, the walker will be encouraged to favor visits of correlated and anti-correlated genes. We have conducted many experiments by replacing the parameters proposed by [START_REF] Mikolov | Distributed representations of words and phrases and their compositionality[END_REF] with our suggested use of similarity between nodes or by combining the different ways to bias the walk. It turns out, in the end, that using only the bias based on the similarity of p-values gives the best results. The bias we introduced allows to control the walk by assigning a transition t from a node i to a node j proportional to:

t i , j =1-| p i -p j |
Other parameters tuning the Node2Vec method are summarized in Table 4.

Algorithm

The cohesion measurement of a set of nodes on the graph can be determined on the embedded space using the cosine distance cos between the vectors representing the nodes. We use this property to identify, the most similar nodes to a given node v i (represented as similar function in the algorithm summarized in Fig. 8). Thanks to a greedy approach, we collect, from each node, clusters M i of increasing size evaluated using the s score previously defined. Our strategy is to expand the cluster as long as the s score increases (lines 8-10 of the algorithm). This allows to overcome the known drawback of the s score which is that high scores obtained by small sub-networks can be overtaken by random scores obtained in large networks [START_REF] Nikolayeva | Network module identification-A widespread theoretical bias and best practices[END_REF]. In practice, as we are very strict on the stopping condition, the clusters obtained are quite small (usually 5 nodes at most). At the end of this phase, we obtain a list of clusters, each one centered on a node, with each cluster being assigned a corresponding s score (Fig. 7D). The cluster centered on vertex v i is thus denoted

M i with M i ⊆ V and v i ∈ M i .
The next phase of the method consists in combining the different clusters while ensuring that new merged clusters remain spatially cohesive (Fig. 7E). We say that two clusters M i and M j are spatially cohesive when M i ∩ M j ≠ ∅ . Starting from the module with the higher score (line 14 of the algorithm), the process consists in evaluating all possible clusters formed by M i ∪ ℘ ( M j / j ≠i ) using the s score, with ℘ ( M i ) denoting the powerset of M i and keeping the modules with the highest s scores (lines 15-20). The workflow of the AMINE method is presented in Fig. 7 and Fig. 8.

Generation of realistic interaction network

We use an extended version of the Barabasi-Albert model of preferential attachment [START_REF] Albert | Topology of evolving networks: local events and universality[END_REF], to generate several artificial networks by varying the parameters p and q controlling the probabilities to add and remove edges respectively as well as the parameter specifying the number of initial nodes. Our results suggest that using 3 initial nodes with parameters p and q set to 0.09 and 0.70 respectively allows to generate random network with topologies relatively close to real interaction networks (Table 5).

1000 graphs were generated using these parameters (an example of this kind of graph is given in supplementary Fig. S2). The strategy to specify the value of nodes is exactly the same at the one of [START_REF] Robinson | Incorporating interaction networks into the determination of functionally related hit genes in genomic experiments with Markov random fields[END_REF]. In order to be able to do a comparison with other methods, we generated only one module of designated hits. As AMINE is dedicated to the identification of relatively small modules (in order to focus on really relevant genes that can be investigated by biologists), we have targeted our tests on the identification of small modules of size 10 and 20.

Cell culture

The PDAC cell line, MIA PaCa-2, was obtained from Richard Tomasini CRCM, Marseille, France and culture in DMEM (Gibco™, Life Technologies Limited, Paiseley, UK) supplemented with 10% FBS, and penicillin/streptomycin. Cells were maintained at 37° C in a humidified atmosphere (5% CO 2 ). Cells were tested routinely for Mycoplasma contamination.

siRNA transfection siRNAs (Sigma-Aldrich) were used for BLIMP 1 silencing. Non-targeting (si-Ctrl: SIC001) or BLIMP 1-targeting siRNAs (si-Blimp 1-1: 5'CUUGGAAGAUCUGACCCGA3'; si-Blimp 1-2: 5'CCUUUCAAAUGUCAGACUU3' were transfected in MIA PaCa-2 cells using Lipofectamine ® RNAiMAX (Invitrogen, Life technologies Corp., Carlsbad, California, USA) following the manufacturer's instructions. The final siRNA concentration was 30 nM. The medium was changed 8 h after transfection and the efficiency of the transfection was assessed by western blot after 72h.

Western blotting

Cells were lysed in RIPA buffer supplemented with Complete Protease Inhibitor Cocktail and PhosSTOP Phosphatase Inhibitor Cocktail (Roche Diagnostics GmbH, Mannheim, Germany).

Lysate were centrifuged at 12 000 rpm for 15 min at 4°C and then protein concentration was quantified using Bradford assay. Protein lysate were subject to SDS-PAGE and transferred onto a PVDF membrane. Membranes were blocked with 5% low fat milk in Tris-buffer saline -tween (TBS-T) for 1 hour. Membrane were incubated in Blimp 1 antibody (Cell Signaling, diluted at 1:1000) overnight. Membrane were washed in TBS-T followed by incubation with Horseradish peroxidase-conjugated secondary antibody for 1 hour at room temperature (Sigma-Aldrich). Signal was then visualized using ECL reagent (Immobilon ® Western, Millipore, Burlington, MA, USA) and chemoluminescence detection system (fusion FX7 Edge; Vilber, Marne-la-Vallée, France).

Human Cytokine Array

For the cytokine assay, the Proteome Profiler™ Human XL Cytokine Array Kit (R&D Systems, Minneapolis, MN, USA) was used. The array was carried out using 500µl of cell supernatants obtained by incubating MIA PaCa-2 in DMEM 0% FBS, 48h after siRNA (si-Ctrl or si-Blimp 1-2) transfection, following the manufacturer's instructions. For analysis of cytokine arrays, the intensity of each spot was measured using ImageJ software. Background was removed from all values, and they were normalized to the positive control spots. Comparison of several metrics associated, on the one hand, an artificial dense network generated with the extended Barabasi-Albert model using 3 initial nodes and setting parameters p and q to 0.09 and 0.7 respectively and, on the other hand, a subnetwork containing the same number of vertices extracted from the STRING database.

Figures and Tables

Supplementary Materials

Figs. S1 to S5.

Fig. 1 .

 1 Fig. 1. Evaluation of AMINE on artificial data generated by Robinson et al. (2017). (A) boxplots of the F1, recall and precision measures. The "trucated recall" displays the recall obtained by setting the maximum size of the 3 best identified modules to 30. (B) repartition of the size of the identified modules.

Fig. 2 .

 2 Fig. 2. Results obtained for the identification of modules of size 10. Comparison of AMINE with COSINE, BioNet and GiGA on a task consisting in identifying a module of size 10 on a dense artificial network with 1000 vertices. The boxplots summarize the results obtained on 1000 different networks. "Expr. value", that is used for baseline, consists in selecting the 10 genes with the lowest pvalues. (A) Boxplots representing the distribution of F1 scores. (B) Boxplots representing the distibution of identified module sizes. (C) Magnification of the distribution of module sizes for GiGA and AMINE.

Fig. 3 .

 3 Fig. 3. Results obtained for the identification of modules of size 20. Comparison of AMINE with COSINE, BioNet and GiGA on a task consisting in identifying a module of size 20 on a dense artificial network with 1000 vertices. The boxplots summarize the results obtained on 1000 different networks. "Expr. value", that is used for baseline, consists in selecting the 20 genes with the lowest pvalues. (A) Boxplots representing the distribution of F1 scores. (B) Boxplots representing the distibution of identified module sizes. (C) Magnification of the distribution of module sizes for GiGA and AMINE.

Fig. 4 .

 4 Fig.4. AMINE reveals modules associated to metastatic process in HMGA2 positive PDAC cells. The network was generated by STRING using the five first modules generated by AMINE from the list of deregulated genes in HMGA2+ PDAC cells available in[START_REF] Chiou | BLIMP1 induces transient metastatic heterogeneity in pancreatic cancer[END_REF]

Fig. 5 .

 5 Fig. 5. BLIMP1 is associated to immune response and inflammation in PDAC cells. (A) BLIMP1-associated module generated by the profiling of genes deregulated in PDAC metastatic HMGA2 positive cell population. (B) Modules 2 and 9 generated by the profiling of genes down-regulated in BLIMP1-silenced PDAC cells. The networks were generated by STRING.

Fig. 6 .

 6 Fig. 6. BLIMP1 silencing modify the cytokine secretion profile in PDAC cells. (A) Immunoblots of Blimp 1 in MIA PaCa-2 cells transfected with a non-targeting siRNA (Si-Ctrl), or with two different siRNA targeting Blimp-1 (Si-Blimp 1-1 and Si-Blimp 1-2). Data are representative of three independent experiments. (B) Soluble cytokine protein expression was assessed using cytokine arrays in si-Ctrl or si-Blimp 1 transfected MIA PaCa-2 (N=2; n = 4). Representative arrays are shown. On the left panel, values from densitometry quantification are shown as a fold change from the control.

Fig. 7 .

 7 Fig. 7. Workflow of the AMINE method. (A) Input data are composed of a table storing the significance of the expression variation of genes between two conditions and a network representing known gene interactions. (B) Data about gene interactions and gene variations are merged to generate an attributed gene network. (C) Nodes belonging to the attributed gene network are mapped to a low-dimensional space through the use of a biased Node2Vec method. (D) Sets of genes that are both cohesive and differentially expressed are identified in the embedded space by maximizing both the s scores of the nodes and the cosine distance cos between the vectors representing the nodes. (E) Redundancy in the content of modules is ruled out by combining sets of nodes obtained in the previous step while ensuring that the result remains spatially cohesive.

Fig. 8 .

 8 Fig. 8. Algorithm of the method AMINE.

Table 1 . Functions associated to the five best modules revealed by AMINE profiling of genes deregulated in PDAC metastatic HMGA2 positive cell population.

 1 Gene enrichment analyses were performed using data from the STRING database.

	Module number	P-value of the module	Function/Signaling pathway	Accession Id	P -value of annotation
	1	2.24e-7	Collagen fibril organisation Extra cellular matrix organisation	GO:0030199 GO:0030198	4.48e-13 4.48e-13
	2	7.47e-6	HIF-1 signaling pathway Central carbon metabolism in cancer	mmu04066 mmu05230	6.43e-10 1.17e-05
			cell adhesion	GO:0007155	5.60e-06
	3	4.8e-5	ECM-receptor interaction Focal adhesion	mmu04512 mmu04510	1.54e-08 2.40e-07
			PI3K-Akt signaling pathway	mmu04151	1.36e-06
	4	6.07e-5	extracellular matrix structural constituent GO:0005201 Extracellular matrix organization MMU-1474244 1.14e-06 1.82e-06
			transmembrane receptor protein tyrosine kinase signaling pathway	GO:0007169	4.04e-16
	5	6.17e-5	regulation of cell migration regulation of peptidyl-tyrosine phosphorylation	GO:0030334 GO:0050730	4.64e-07 5.60e-07
			angiogenesis	GO:0001525
	Module number	P-value of the module	Function/Signaling pathway	Accession Id	P -value of annotation
			regulation of interleukin-2 production	GO:0032663	0.0129
	169	0.044	regulation of alpha-beta T cell differentiation	GO:0046637	0.0129
			Th17 cell differentiation	mmu04659	0.00038

Table 2 . Functions associated to the module embedding BLIMP1 revealed by AMINE profiling of genes deregulated in PDAC metastatic HMGA2 positive cell population.

 2 Gene enrichment analyses were performed using data from the STRING database.

	Module number	Function / Signaling pathway	Accession Id	P-value of annotation
		Inflammatory response	GO:0006954	1.05e-12
	2 + 9	cytokine activity	GO:0005125	2.45e-12
		immune system process	GO:0002376	1.73e-12

Table 3 . Functions associated to modules 2 and 9 resulting from AMINE profiling of genes down-regulated in PDAC cells silenced for Blimp1.

 3 The two modules were merged and subjected to gene enrichment analyses using data from the STRING database.

	Parameters	Values
	Number of walks	20
	Walk length	100
	Vector dimensions	128
	Window size	5
	epoch	10

Table 4 . Parameters used for the Node2Vec method.

 4

	Metric	Artificial dense network	Real network
	Number of vertices	5980	5980
	Alpha coefficient	1.908	1.602
	R square	0.876	0.924
	Mean of neighbor connectivity 42.207	47.720

Table 5 . Comparison of artificial and real biological network.

 5

a clique is a subset of vertices of an undirected graph such that every two distinct vertices in the clique are connected[START_REF] Seidman | A graph-theoretic generalization of the clique concept[END_REF].

Acknowledgments

The authors are grateful to the OPAL infrastructure from Université Côte d'Azur and the Université Côte d'Azur's Center for High-Performance Computing for providing resources and support.

Funding UCAJEDI Investments in the Future project managed by the National Research Agency (ANR) under reference number ANR-15-IDEX-01 French National Research Agency (ANR) through the LABEX SIGNALIFE program (reference # ANR-11-LABX-0028-01).

Data and materials availability

All data are available in the main text or the supplementary materials.

Author contributions

Competing interests

Authors declare that they have no competing interests.