Brian Kinney

Shubham Sahu

Natalia Stec

Kelly Hills-Muckey

Dexter W Adams

Matt Jaremko, Leemor Jing Wang

Wolfgang Keil
email: wolfgang.keil@curie.fr

Christopher M Hammell
email: chammell@cshl.edu

Matt Jaremko

Joshua Leemor

-Tor

A circadian-like gene network programs the timing and dosage of heterochronic miRNA transcription during C. elegans development

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

Introduction

Many gene regulatory networks encode and read information controlling the tempo and order of developmental events. C. elegans postembryonic maturation is compartmentalized into four larval stages with distinct patterns of cell division, cell differentiation, and cuticle formation that are separated by molts [START_REF] Sulston | Post-embryonic cell lineages of the nematode, Caenorhabditis elegans[END_REF] . The transition from one stage-specific pattern of cell division to the next is mediated by the accumulation of heterochronic microRNAs (miRNAs) that post-transcriptionally downregulate temporal identity target genes that define stage-specific gene expression patterns [START_REF]MicroRNAs and developmental timing[END_REF] .

Hence, lower or higher dosages of heterochronic miRNAs result in precocious cellular transitions or reiteration of stage-specific cellular programs, respectively. The transcription of heterochronic microRNAs is periodic, peaking at specific phases of each larval molting cycle [START_REF] Van Wynsberghe | Period homolog LIN-42 regulates miRNA transcription to impact developmental timing[END_REF][START_REF] Perales | the Caenorhabditis elegans PERIOD homolog, Negatively Regulates MicroRNA Transcription[END_REF][START_REF] Kim | Dampening of expression oscillations by synchronous regulation of a microRNA and its target[END_REF][START_REF] Mcculloch | Caenorhabditis elegans period homolog lin-42 regulates the timing of heterochronic miRNA expression[END_REF] . This suggests the existence of a clock-like system to coordinate their temporal expression with overall developmental progression.

Evidence suggests that LIN-42, the C. elegans ortholog of the circadian Period protein, is required for temporal patterning by directly modulating the dynamic features of heterochronic microRNA transcription [START_REF] Van Wynsberghe | Period homolog LIN-42 regulates miRNA transcription to impact developmental timing[END_REF][START_REF] Perales | the Caenorhabditis elegans PERIOD homolog, Negatively Regulates MicroRNA Transcription[END_REF][START_REF] Mcculloch | Caenorhabditis elegans period homolog lin-42 regulates the timing of heterochronic miRNA expression[END_REF] . Loss-of-function mutations in lin-42 increase the amplitude and duration of cyclical heterochronic miRNA transcription, indicating that LIN-42, like its mammalian ortholog, functions as a transcriptional repressor [START_REF] Van Wynsberghe | Period homolog LIN-42 regulates miRNA transcription to impact developmental timing[END_REF][START_REF] Perales | the Caenorhabditis elegans PERIOD homolog, Negatively Regulates MicroRNA Transcription[END_REF][START_REF] Mcculloch | Caenorhabditis elegans period homolog lin-42 regulates the timing of heterochronic miRNA expression[END_REF] . Consequently, lin-42 mutants develop precociously due to overaccumulation of heterochronic microRNAs [START_REF] Perales | the Caenorhabditis elegans PERIOD homolog, Negatively Regulates MicroRNA Transcription[END_REF][START_REF] Jeon | Similarity of the C. elegans developmental timing protein LIN-42 to circadian rhythm proteins[END_REF] . While most core components of the heterochronic gene regulatory network are expressed in a graded temporal pattern, LIN-42 is highly dynamic, with a single peak of expression during each larval stage [START_REF] Perales | the Caenorhabditis elegans PERIOD homolog, Negatively Regulates MicroRNA Transcription[END_REF][START_REF] Jeon | Similarity of the C. elegans developmental timing protein LIN-42 to circadian rhythm proteins[END_REF] . The similar molecular functions and expression dynamics between LIN-42 and mammalian Period suggest that a regulatory architecture akin to the one that generates the once-daily circadian transcriptional patterns in mammals may play a role in temporal cell fate specification in nematodes. Intriguingly, C. elegans lacks orthologs of CLOCK and BMAL1, the central transcription factors that drive oscillatory circadian transcription and the direct targets of Period repression. Thus, the mechanism by which LIN-42 modulates the transcriptional output of C. elegans heterochronic microRNAs is currently unknown.

Oscillatory lin-4 transcription is pulsatile

Previous measurements of miRNA transcription used destabilized GFP reporters that must be transcribed, processed, and translated to visualize expression dynamics [START_REF] Perales | the Caenorhabditis elegans PERIOD homolog, Negatively Regulates MicroRNA Transcription[END_REF][START_REF] Stec | An Epigenetic Priming Mechanism Mediated by Nutrient Sensing Regulates Transcriptional Output during C. elegans Development[END_REF] . These features limit their temporal resolution. We used the MS2/MCP-GFP tethering system, where engineered RNA loops derived from MS2 bacteriophage and a co-expressed MS2 coat protein fused to green 4 fluorescent protein (MCP-GFP) can be concentrated and localized at a gene of interest by active transcription (Fig. 1a) [START_REF] Tutucci | An improved MS2 system for accurate reporting of the mRNA life cycle[END_REF] . We measured the transcriptional dynamics of the lin-4 heterochronic miRNA that is expressed periodically throughout all larval stages and down-regulates lin-14 expression within the heterochronic gene regulatory network early in development [START_REF] Stec | An Epigenetic Priming Mechanism Mediated by Nutrient Sensing Regulates Transcriptional Output during C. elegans Development[END_REF] . We generated a transgene harboring 24 copies of a synthetic MS2 hairpin immediately downstream of the lin-4 pre-miRNA (Fig. 1a and Extended Data Fig. 1a) [START_REF] Stec | An Epigenetic Priming Mechanism Mediated by Nutrient Sensing Regulates Transcriptional Output during C. elegans Development[END_REF][START_REF] Bracht | Regulation of lin-4 miRNA expression, organismal growth and development by a conserved RNA binding protein in C. elegans[END_REF] . This transgene rescues lin-4 null allele phenotypes (Extended Data Fig. 1b,c). We also ubiquitously expressed MCP-GFP to detect MS2-tagged RNAs and a histone::mCherry fusion to locate nuclei (Extended Data Fig. 1a). Examination of transgenic animals revealed that the formation of nuclear MCP-GFP foci occurred at each larval stage in somatic tissue types known to transcribe lin-4 (Fig. 1b). MCP-GFP foci were not observed in developing embryos (n > 50) or in starvation arrested L1 larva (n = 23); consistent with the activation of lin-4 transcription after the initiation of larval development 14 .

To directly determine the level of coordination of transcriptional programs, at single-cell resolution, across cell and tissue types in living animals, we used our microfluidics-based platform [START_REF] Keil | Long-Term High-Resolution Imaging of Developing C. elegans Larvae with Microfluidics[END_REF] for long-term imaging of C. elegans larvae harboring lin-4::24xMS2/MCP-GFP system. We quantified expression dynamics in hypodermal cells where the timing of transcriptional activity can be accurately assessed in relation to the occurrence of stage-specific cell division patterns [START_REF] Sulston | Post-embryonic cell lineages of the nematode, Caenorhabditis elegans[END_REF] . We screened for periods of lin-4 transcriptional activity by imaging at 15min time intervals from the first larval stage (L1) to the mid-L4 stage (~60h) (n>10). lin-4::24xMS2 transcription was highly pulsatile, with a single transcriptional episode of ~45-90min in each cell during each larval stage, followed by long periods of inactivity (>8hrs at 20°C) (Fig. 1c). Transcriptional activation across cells within the hypodermis was highly concordant, exhibiting similar transcriptional on and off times for all lin-4::24xMS2 loci (Fig. 1e). For instance, during the mid-L3 stage, the appearance of lin-4::24xMS2 expression throughout the hypodermis generally occurred within minutes of the first vulval precursor cell (VPC) divisions (P3.p or P4.p) and were completed by the divisions of remaining VPCs (P5.p-P7.p) (n = 15) (Fig. 1c).

To determine how similar the transcriptional dynamics in individual cells within the hypodermis are, we performed short-term imaging time courses (<6h) at 4min intervals in staged larvae. During transcriptional episodes, we detected near synchronous accumulation of MCP-GFP foci at each hypodermal lin-4::24xMS2 locus for 60-90 minutes (Fig. 1d,e, Suppl. Movie 1,2) (>15 animals). We found no signs of "bursty" transcription [START_REF] Larson | A single molecule view of gene expression[END_REF] as MCP-GFP foci were continuously maintained at each lin-4::24xMS2 transgene for the entire transcriptional episode (Suppl. Movie 1,2). These features were independent of the number of lin-4::24xMS2 loci per nucleus as cell types that undergo endoreduplication (i.e., hyp7 cells) exhibited MCP-GFP foci dynamics indistinguishable from diploid cells (Fig. 1d). The dynamic features of lin-4::24xMS2 expression in hypodermal cells were similar across different developmental stages (Fig. 1e) suggesting that the same regulatory programs controlling lin-4 transcription were repeated at each larval stage. Pulsatile transcription also occurs in a broad array of additional cell types that normally express lin-4 (including the pharynx and intestinal cells) (Extended Data Fig. 2). Therefore, the gene regulatory network that generates lin-4 transcriptional pulses at each stage of development organizes the timing, amplitude and duration of transcription in a highly reproducible manner.

NHR-85 and NHR-23 function cooperatively

Full transcriptional activation of lin-4 requires a conserved upstream regulatory element, the Pulse Control Element (PCE), located ~2.8kb upstream of the lin-4 sequence 11 . We performed a yeast-one-hybrid screen to discover TFs that bind this regulatory region using the entire 514bp PCE as bait. We identified three TFs that specifically bound the PCE (BLMP-1, NHR-23, and NHR-85) (Fig. 2a). BLMP-1 controls the amplitude of lin-4 expression and functions as a pioneer factor to decompact the lin-4 locus throughout development [START_REF] Stec | An Epigenetic Priming Mechanism Mediated by Nutrient Sensing Regulates Transcriptional Output during C. elegans Development[END_REF] . nhr-85 and nhr-23 encode two nuclear hormone receptors (NHRs) that are the closest nematode orthologs of human circadian TFs Rev-Erb and ROR, respectively (Fig. 2b). NHR-23 expression cycles with the larval molts and has been implicated in controlling the pace of larval development [START_REF] Kostrouchova | Nuclear hormone receptor CHR3 is a critical regulator of all four larval molts of the nematode Caenorhabditis elegans[END_REF][START_REF]d) Time course experiments demonstrate that lin-4::24xMS2 expression occurs immediately before the first Pn[END_REF][START_REF] Macneil | Diet-induced developmental acceleration independent of TOR and insulin in C. elegans[END_REF] . Furthermore, nhr-23 and the terminal heterochronic miRNA, let-7, genetically interact to mutually limit supernumerary molts during adulthood [START_REF]d) Time course experiments demonstrate that lin-4::24xMS2 expression occurs immediately before the first Pn[END_REF][START_REF] Hayes | The mir-84 and let-7 paralogous microRNA genes of Caenorhabditis elegans direct the cessation of molting via the conserved nuclear hormone receptors NHR-23 and NHR-25[END_REF] . Analysis of publicly available ChIP-seq data indicated that all three TFs interact in vivo with lin-4 regulatory sequences (Fig. 2a), and their binding sites are enriched in the promoters of cyclically expressed mRNAs and other heterochronic miRNAs (Extended Data Fig. 3) (table S1 andS2).

Nuclear hormone receptors often bind cooperatively as homo-or hetero-dimeric complexes at closely spaced cis-regulatory DNA elements [START_REF] Forman | Unique response pathways are established by allosteric interactions among nuclear hormone receptors[END_REF][START_REF] Perlmann | Two distinct dimerization interfaces differentially modulate target gene specificity of nuclear hormone receptors[END_REF][START_REF] Rastinejad | Structural determinants of nuclear receptor assembly on DNA direct repeats[END_REF] . Several features of NHR-85 Rev-Erb , NHR-23 ROR, and the lin-4 PCE suggest that this may also be the case for lin-4 transcription. First, NHR-85 Rev-Erb and NHR-23 ROR share significant sequence homology within their C4-type Zinc finger DNA-binding domains, suggesting they bind similar sites (Fig. 2c) [START_REF] Burris | Nuclear hormone receptors for heme: REV-ERBalpha and REV-ERBbeta are ligand-regulated components of the mammalian clock[END_REF] . Second, we found that NHR-85 Rev-Erb and NHR-23 ROR heterodimerize in yeast two-hybrid assays and in vitro with high affinity (5.8 +/-2.2 nM KD (Fig. 2d,e)). Third, sequences within the PCE element contain direct GGTCA repeats that are predicted binding sites of NHR-85 Rev-Erb and NHR-23 ROR (Fig. 2f andg). Fourth, we found that concentrations of NHR-23 ROR that were insufficient to bind the PCE alone were dramatically stimulated by the addition of NHR-85 Rev-Erb , indicating cooperative binding (Fig. 2h,i). Supporting this, 6 NHR-85 Rev-Erb alone could not bind lin-4 PCE DNA fragments (Extended Data Fig. 4c). Finally, the cooperative binding of NHR-85 Rev-Erb and NHR-23 ROR required both GGTCA repeats, suggesting they bind closely spaced regulatory elements as a heterodimer (Fig. 2h).

NHR-85/NHR-23 expression coincides with lin-4 transcription

We quantified mRNA and protein expression during post-embryonic development to determine how NHR-85, NHR-23, and LIN-42 expression patterns may contribute to the regulation of oscillatory lin-4 transcription. Expression of nhr-85 begins from an L1-stage arrest with a pulse of transcription that is followed by a monotonic expression pattern for the remaining larval stages (Fig. 3a). In contrast, nhr-23 and lin-42 mRNAs were expressed in phased, high-amplitude oscillatory patterns (Fig. 3a). We next explored the temporal dynamics of the corresponding proteins, by quantifying the expression of endogenously-tagged alleles during the L4 stage, where changes in vulval morphogenesis can be directly correlated with developmental age [START_REF] Mok | Morphologically defined sub-stages of C. elegans vulval development in the fourth larval stage[END_REF] . NHR-23 ROR ::mScarlet and LIN-42 Period ::GFP are dynamically expressed in all hypodermal cells, with a single peak of expression that matched the phased expression of their mRNAs (Fig. 3c). We also found that NHR-85 Rev-Erb ::GFP expression was highly dynamic during these periods, indicating substantial posttranscriptional regulation of expression in the L2-L4 stages of development. Specifically, peak expression of NHR-85 Rev-Erb ::GFP occurred at ecdysis (shortly before NHR-23 ROR ::mScarlet onset), becomes undetectable by the L4.3 stage of development, and resumes expression at the L4.6 stage in an antiphasic manner to the expression pattern of LIN-42 Period ::YFP (Fig. 3b). Highly similar phased expression patterns of these proteins are also maintained in L4-staged vulval tissues (Fig. 3c) as well as in L3-staged hypodermal cells and VPCs (Fig. 3c-e).

Since NHR-85 Rev-Erb and NHR-23 ROR heterodimerize and bind cooperatively to regions of the lin-4 enhancer that control dynamic transcription, we hypothesized that the 60-90min pulses of lin-4 transcription might only occur in the short window of each larval stage where NHR-85 Rev-Erb and NHR-23 ROR are co-expressed. To compare the timing of these events, we examined MCP-GFP localization in vulval cells during the L3 and L4 stages, where the rapid, stereotyped vulval cell division patterns [START_REF] Keil | Long-Term High-Resolution Imaging of Developing C. elegans Larvae with Microfluidics[END_REF] and changes in morphology [START_REF] Mok | Morphologically defined sub-stages of C. elegans vulval development in the fourth larval stage[END_REF] enable precise determination of the timing of lin-4::24xMS2 transcription and TF expression dynamics. We found a correspondence between NHR-85 Rev-Erb and NHR-23 ROR co-expression and lin-4::24xMS2 transcription in both vulval and hypodermal cells (Fig. 3d,e). While the expression of both nuclear receptors is phased, the transient expression of the lin-4::24xMS2 transgene only occurs during the brief period when both NHRs are expressed. The timing of NHR-85 Rev-Erb post-transcriptional downregulation is also concurrent with the onset of LIN-42 Period expression in VPCs and hypodermal cells (Fig. 3c-e). This suggests that the dynamic patterns of these TFs control the timing and duration of lin-4 transcriptional pulses.

LIN-42 Period posttranscriptionally represses NHR-85

Human Per2 (PERIOD2) protein interacts with multiple mammalian NHRs (including Rev-Erb) to modulate their transcriptional activity [START_REF] Schmutz | The mammalian clock component PERIOD2 coordinates circadian output by interaction with nuclear receptors[END_REF] . To test whether LIN-42 Period physically interacts with either NHR-85 Rev-Erb or NHR-23 ROR , we used two-hybrid assays [START_REF] Walhout | Yeast two-hybrid systems and protein interaction mapping projects for yeast and worm[END_REF] . We found that both major LIN-42 isoforms interact with NHR-85 Rev-Erb but not NHR-23 ROR (Fig. 4a). We mapped the regions of LIN-42 Period that are required for NHR-85 Rev-Erb binding and found that a minimal 51aa fragment present in both major LIN-42 Period isoforms is sufficient to mediate interactions (Fig. 5a). This domain is distinct from the interaction motifs implicated in mammalian Per2 and Rev-Erb 26 . Because Per2 has been shown to interact with many mammalian NHRs, we extended our two-hybrid analysis of LIN-42 Period interactors to determine if LIN-42 Period also binds additional NHRs. We performed two-hybrid experiments between LIN-42 Period isoforms and 241 of the remaining 282 encoded C. elegans NHRs. We identified 65 NHRs that physically interact with LIN-42 Period (Extended Data Fig. 5a,b). The additional interacting NHRs included DAF-12, which regulates the expression of the let-7-family of miRNAs and controls dauer development [START_REF] Tennessen | elegans developmental timing protein LIN-42 regulates diapause in response to environmental cues[END_REF][START_REF] Hammell | A feedback circuit involving let-7-family miRNAs and DAF-12 integrates environmental signals and developmental timing in Caenorhabditis elegans[END_REF][START_REF] Bethke | Nuclear hormone receptor regulation of microRNAs controls developmental progression[END_REF] , and NHR-14 HNF4a , NHR-69 HNF4a , and NHR-119 PPARa (Extended Data Fig. 5a,b) whose orthologs are also bound by Per2 [START_REF] Schmutz | The mammalian clock component PERIOD2 coordinates circadian output by interaction with nuclear receptors[END_REF] . These findings suggest that many physical interactions between Period orthologs have been maintained since the divergence of nematodes and man and are, therefore, likely functional.

Given the physical interaction between NHR-85 Rev-Erb and LIN-42 Period , we asked whether NHR-85 Rev-Erb expression was required for the precocious phenotypes seen in lin-42(lf) mutants (lin-42(n1089)). We found that lin-42(lf) heterochronic phenotypes are partially ameliorated by removing nhr-85 function. Specifically, the precocious expression of adult-specific reporters (e.g., col-19::GFP) in both seam and hyp7 cells observed in lin-42(lf) mutants is suppressed by nhr-85 deletion; leaving weak expression in seam cells in double mutants, while precocious deposition of adult alae was not suppressed (Fig. 4b). To examine whether LIN-42 Period regulates NHR-85 Rev-Erb temporal expression, we compared the dynamics of NHR-85 Rev-Erb ::GFP (and NHR-23 ROR ::mScarlet) levels in wild-type and lin-42(lf) mutants. We found that the expression of NHR-85 Rev-Erb is altered in two ways by the lin-42 mutation. First, the expression of NHR-85 Rev-Erb ::GFP is ~2.3x more abundant at the beginning of the L4 stage in lin-42 mutants when compared to wild-type animals (Fig. 4c,d). More importantly, the periodic dampening of NHR-85 Rev-Erb expression that usually occurs by the L4.2 stage of vulval morphogenesis (Fig. 3c,d) is altered in lin-42(lf) mutants. Specifically, NHR-85 Rev-Erb ::GFP expression perdures into the L4.4 stage lin-42(lf) animals (Fig. 4c,d). Mutations in lin-42 do not alter the onset or duration of NHR-23 ROR ::mScarlet accumulation in hypodermal or vulval cells. This suggests that LIN-42 regulates lin-4 transcriptional output by controlling the duration of NHR-85 Rev-Erb /NHR-23 ROR heterodimeric complex formation in a manner that directly correlates with NHR-85 abundance.

Given the above hypothesis, we anticipated that specific features of lin-4 transcriptional pulses would be altered in lin-42 mutants. While most wild-type seam cells exhibit MCP-GFP foci during each larval stage (76%; n=66 nuclei), the percentage of seam cells showing detectable lin-4 transcription is dramatically increased in lin-42 mutants (100%; n=58 nuclei). This indicates that lin-42 normally dampens lin-4 transcriptional pulses in wild-type animals. In addition to elevating the likelihood that lin-4::24xMS2 transcription is above a threshold sufficient to generate measurable MCP-GFP foci, lin-42(lf) mutations lead to an increase in the intensity of MCP-GFP foci indicating that LIN-42 Period normally also limits the rate of transcriptional activation of the lin-4 locus (Fig. 4e).

Time course experiments also revealed that overall duration of transcription events in lateral seam cells was ~2.2 times longer in lin-42 mutants compared to wild-type (Fig. 4e). Consistent with our hypothesis that the duration of lin-4 transcription is dependent on the duration of NHR-85 Rev-Erb expression, the transcriptional onset of lin-4::24xMS2 occurs earlier in L3 staged lin-42 mutants and inappropriately extended through both the first and second Pn.p divisions (Fig. 4f).

While nhr-85(0) mutants did not misexpress adult-specific reporters (Fig. 5b), high-resolution imaging of VPC divisions and lin-4::24xMS2 expression indicate two features of developmental timing were altered. In wild-type animals, transcriptional pulses of lin-4::24xMS2 are both robust and concordant in adjacent wild-type VPCs (Figs. 4g &1c). In contrast, MCP-GFP foci in nhr-85(0) mutants begin to accumulate at the same relative phase of L3-stage VPC development but are dimmer and only transiently observed (Fig. 4g). Second, under identical imaging conditions, the rapid and highly coordinated VPC divisions observed in wild-type animals is altered in nhr-85(0) mutants with some P5.p and P7.p dividing hours after the first P6.p division (Fig. 4g). These results indicate that NHR-85 functions to enhance the robustness of temporally regulated processes during development and that some level of lin-4 transcription occurs without NHR-85 Rev-Erb , perhaps driven by NHR-23 ROR alone. RNAi-mediated depletion of nhr-23 activity in wild-type animals resulted in mild heterochronic phenotypes (Extended Data Fig. 6). Consistent with the hypothesis that NHR-23 and NHR-85 function cooperatively to control temporal regulation, the penetrance of these phenotypes was enhanced when nhr-23 was also depleted in nhr-85(0) animals (Extended Data Fig. 6).

Discussion

This study reveals the gene regulatory network controlling lin-4 transcriptional pulses. This network shares integral components with the human circadian clock but exhibits essential differences in its regulatory architecture. In the circadian clock regulatory network, transcription factors CLOCK and BMAL1 generate rhythmic expression patterns of clock control genes (CCGs) [START_REF] Bunger | Mop3 is an essential component of the master circadian pacemaker in mammals[END_REF][START_REF] Hogenesch | The basic-helix-loop-helix-PAS orphan MOP3 forms transcriptionally active complexes with circadian and hypoxia factors[END_REF][START_REF] King | Positional cloning of the mouse circadian clock gene[END_REF][START_REF] Reick | NPAS2: an analog of clock operative in the mammalian forebrain[END_REF] (Fig. 5h), including two core transcriptional repressors, Period and CRY as well as NHR genes Rev-Erb and ROR. Negative feedback on CLOCK/BMAL1 expression by Period/CRY heterodimers is essential for the proper generation of circadian rhythms. Rev-Erb and ROR modulate CLOCK and BMAL1 expression through opposing transcriptional activities (Fig. 4h) [START_REF] Preitner | The orphan nuclear receptor REV-ERBalpha controls circadian transcription within the positive limb of the mammalian circadian oscillator[END_REF][START_REF] Ueda | A transcription factor response element for gene expression during circadian night[END_REF][START_REF] Sato | A functional genomics strategy reveals Rora as a component of the mammalian circadian clock[END_REF] but are dispensable for the generation of circadian oscillations 35 . Instead of CLOCK and BMAL1, which are lacking in the C.

elegans genome, we propose here that the worm orthologs of Rev-Erb and ROR compose the central transcription factors of the hypodermal developmental clock (Fig. 5h). In contrast to their antagonistic roles in the circadian clock, C. elegans NHR-85 Rev-Erb and NHR-23 ROR heterodimerize and work cooperatively to induce transcriptional pulses (Fig. 4h). The evolutionary rewiring of this interaction results in the phased expression of NHR-85 Rev-Erb and NHR-23 ROR generating a short temporal window for these proteins to cooperativity promote transcription. This enables transcription to be precisely dosed and timed at each larval stage. This rewiring also provides a direct physical link to the conserved opposing arm of the clock, LIN-42 Period , which binds to NHR-85 Rev-Erb to control NHR-85 Rev- Erb expression post-transcriptionally (Fig. 4h). Additional physical and functional interactions between Period orthologs and conserved NHRs are maintained in the mammalian and C. elegans systems [START_REF] Schmutz | The mammalian clock component PERIOD2 coordinates circadian output by interaction with nuclear receptors[END_REF][START_REF] Tennessen | elegans developmental timing protein LIN-42 regulates diapause in response to environmental cues[END_REF] and suggest that C. elegans NHRs, in addition to NHR-85 Rev-Erb , may also function in temporal patterning.

While NHR-85

Fig. 1 .

 1 Fig. 1. lin-4 transcription is highly pulsatile at each larval stage. (a) The MS2/MCP system is composed of an MS2 coat protein GFP fusion (MCP-GFP) which can bind to MS2 RNA hairpins engineered into primary miRNA transcripts. (b) Image of mid-L3 staged larva and magnified insets showing MCP-GFP spots in hypodermal cells along the anteroposterior axis. (c) Two examples of L3/L4 hypodermal and vulval precursor cell lineages, with overlayed expression patterns indicating when MCP-GFP foci were visible in each lineage. Red boxes show the temporal region where coordinated lin-4 transcription is observed. Dashed horizontal grey bars indicate hours. (d) Snapshots of individual seam cell and hyp7 cell expression trajectories from L3-staged animals. (e) Expression traces in pairs of V6.p seam cells from L2, L3, and L4 staged animals.

Fig. 2 .

 2 Fig. 2. NHR-85 Rev-Erb and NHR-23 ROR form a heterodimeric complex that binds cooperatively to regulatory elements controlling pulsatile lin-4 transcription. (a) BLMP-1, NHR-85, and NHR-23 bind to the lin-4 PCE in one-hybrid assays 11 . Browser tracks showing BLMP-1, NHR-23, and NHR-85 bindings sites near the lin-4 locus. Also indicated are the major lin-4 pri-miRNAs 13 , the computationally defined binding sites for each TF, and the sub-fragments of the lin-4 PCE element used in gel shifts below. The asterisk indicates the location of the direct repeats of GGTCA in the PCE element. (b) Sequence relationships between NHR-23 and NHR-85 and human Rev-Erb and ROR. (c) NHR-23 ROR and NHR-85 Rev-Erb domain organization. (d) NHR-23 ROR and NHR-85 Rev-Erb interact with each other in two-hybrid assays. (e) Recombinant, purified strep-NHR-85 Rev-Erb, and NHR-23 ROR interact in MST binding assays. (f) Sequences of the PCE sub-fragment iii from a and the 51bp minimal binding element (grey box) derived from lin-4 PCEiii. (g) Conservation of the direct repeats DRs in the PCE elements of different nematode species. (h and i) EMSA experiments of wildtype and mutant target DNAs using recombinant NHR-85 Rev-Erb and NHR-23 ROR .

Fig. 3 .

 3 Fig. 3. lin-4 transcription occurs during the brief periods where NHR-85 Rev-Erb ::GFP and NHR-23 ROR ::mScarlet are co-expressed. (a) RNA-seq time course data of nhr-85, nhr-23, and lin-42 mRNA expression patterns 3 . (b-c) Quantification and micrographs depicting NHR-85 Rev-Erb ::GFP, NHR-23 ROR ::mScarlet, and LIN-42 Period ::YFP expression in hypodermal seam cells and vulval cells, respectively, in each morphologically defined L4 substage 25 . Circles in b represent average

Fig. 4 .

 4 Fig. 4. LIN-42 Period binds to and regulates the expression dynamics of NHR-85 to control the amplitude and duration of lin-4 transcription. (a) LIN-42 isoforms interact with NHR-85 Rev-Erb but not NHR-23 ROR in two-hybrid assays. (b) lin-42(lf) mutants express col-19::GFP during the L3 stage of development and deletion of nhr-85(0) suppresses these phenotypes. Yellow arrows indicate the

  Rev-Erb /NHR-23 ROR heterodimeric complexes generate pulses of lin-4 transcription, NHR-23 ROR may play a more pervasive role in coordinating periodic transcription with overall animal development. Unlike nhr-85, nhr-23 is an essential gene, and nhr-23(0) worms arrest during late embryogenesis/hatching, a time when oscillatory transcription begins[START_REF] Kostrouchova | Nuclear hormone receptor CHR3 is a critical regulator of all four larval molts of the nematode Caenorhabditis elegans[END_REF][START_REF] Meeuse | Developmental function and state transitions of a gene expression oscillator in Caenorhabditis elegans[END_REF][START_REF] Hendriks | Extensive Oscillatory Gene Expression during C. elegans Larval Development[END_REF] . Depletion of NHR-23 ROR during post-embryonic development results in highly penetrant larval arrest phenotypes[START_REF] Johnson | NHR-23 activity is necessary for developmental progression and apical extracellular matrix structure and function[END_REF] . These arrests prevent somatic cell proliferation/differentiation, occur at the beginning of each larval stage, and correlate with the developmental periods where NHR-23 ROR expression peaks[START_REF] Johnson | NHR-23 activity is necessary for developmental progression and apical extracellular matrix structure and function[END_REF] . While the developmental arrest in NHR-23ROR -depleted animals makes it unfeasible to measure changes in transcription for individual genes, these arrests resemble post-embryonic developmental checkpoints where somatic cell proliferation and cyclical gene expression patterns (including lin-4) are halted during acute food removal/starvation[START_REF] Stec | An Epigenetic Priming Mechanism Mediated by Nutrient Sensing Regulates Transcriptional Output during C. elegans Development[END_REF][START_REF] Schindler | Identification of late larval stage developmental checkpoints in Caenorhabditis elegans regulated by insulin/IGF and steroid hormone signaling pathways[END_REF] . Starvation-induced checkpoints are mediated through the 10 regulation of the conserved insulin signaling pathway, which generates various sterol-derived hormones[START_REF] Baugh | Starvation Responses Throughout the Caenorhabditis elegans Life Cycle[END_REF] . These hormones may serve as ligands for individual NHRs and provide a mechanism to organize global gene expression patterns. This type of coordination may underly one mechanism C. elegans larvae employ to generate identical developmental outcomes in diverse and rapidly changing environments.


Maciej Bednarczyk

Hassan Omran

Bernard Bayle

Model Predictive Impedance Control

Keywords: Impedance control, collaborative robotics, physical human-robot interaction

Robots are more and more often designed in order to perform tasks in synergy with human operators. In this context, a current research focus for collaborative robotics lies in the design of high-performance control solutions, which ensure security in spite of unmodeled external forces. The present work provides a method based on Model Predictive Control (MPC) to allow compliant behavior when interacting with an environment, while respecting practical robotic constraints. The study shows in particular how to define the impedance control problem as a MPC problem. The approach is validated with an experimental setup including a collaborative robot. The obtained results emphasize the ability of this control strategy to solve constraints like speed, energy or jerk limits, which have a direct impact on the operator's security during human-robot compliant interactions.

I. INTRODUCTION

The robotics community has been developing for many years now the idea that novel robotic systems should be designed and controlled to perform tasks in contact with their environment, and in particular with human operators. This trend has been increased by the development of collaborative robots. This new type of robot is designed to increase the physical capabilities of human operators, or limit their fatigue. A subsequent research effort has already been devoted to introduce collaborative robots in factories [START_REF] Goodrich | Human-robot interaction: A survey[END_REF], [START_REF] Ajoudani | Progress and prospects of the human-robot collaboration[END_REF]. One of the main challenges is to provide control laws allowing convenient and safe interactions between the human, the robot and the environment [START_REF] Buerger | Complementary stability and loop shaping for improved human-robot interaction[END_REF].

In robotics, when no interaction with the environment is needed, motion control strategies are used. At the contrary, in the presence of interactions, direct force control strategies are preferred for fine force tracking. However, they require good models of the interaction and of the environment, and are not compatible with unpredictable interactions [START_REF] Siciliano | Springer Handbook of Robotics[END_REF]. For this reason, interaction control offers a compromise in order to deal with both force and motion during interactions. Impedance control (IC) [START_REF] Hogan | Impedance control: An approach to manipulation[END_REF] is a widely used interaction control method, particularly efficient for human-robot interactions. It consists in imposing an impedance model for the relationship between the manipulator and the environment. Since its introduction by Hogan, IC has attracted researchers' attention, leading to improvements in handling robots flexibility [START_REF] Ott | Cartesian impedance control of redundant and flexible-joint robots[END_REF] or variable parameters [START_REF] Kronander | Stability considerations for variable impedance control[END_REF], to mention just a few.

In addition to the interaction management, the robot controller should handle several constraints, such as speed and energy limits [START_REF] Lasota | A survey of methods for safe humanrobot interaction[END_REF], limited jerk [START_REF] Macfarlane | Jerk-bounded manipulator trajectory planning: design for real-time applications[END_REF] or actuators saturation. 1 ICube, UMR 7357, University of Strasbourg, CNRS, 1 place de l'Hôpital, 67091, Strasbourg, France. m.bednarczyk@unistra.fr These limits allow ensuring safety [START_REF] Haddadin | Physical Safety in Robotics[END_REF] as well as good ergonomics of the human-robot interaction [START_REF] Kim | Anticipatory robot assistance for the prevention of human static joint overloading in human-robot collaboration[END_REF]. However, it is not straightforward to respect these constraints using classical IC. Model Predictive Control (MPC) is an advanced control method to deal with multi-variable systems with multiple constraints. First proposed for process control in the petrochemical industry [START_REF] Richalet | Algorithmic control of industrial processes[END_REF], MPC has successfully replaced classical methods in many applications [START_REF] Qin | A survey of industrial model predictive control technology[END_REF]. By measuring the current state of the system, MPC predicts the evolution of the system variables over some horizon and computes the optimal input with respect to some performance index, while respecting imposed constraints.

MPC has been only recently used for interaction control. In [START_REF] Killpack | Model predictive control for fast reaching in clutter[END_REF], MPC is used as a high-level controller for path set-point generation, taking into account contact forces due to interactions with the environment. The authors of [START_REF] Matschek | Force feedback and path following using predictive control: Concept and application to a lightweight robot[END_REF] investigate a direct force control strategy based on predictive path following using MPC. In [START_REF] Wahrburg | MPC-based admittance control for robotic manipulators[END_REF], MPC is combined with admittance control for dealing with stability issues while interacting with a stiff environment. Even more recently, nonlinear MPC is formulated by adding admittance dynamics into a path following problem [START_REF] Kazim | Combined predictive path following and admittance control[END_REF]. However, in all these applications, the choice of the cost functions is not directly related to the desired compliant behavior.

To the best of our knowledge, no research has yet been done to combine IC and MPC advantages. Our contribution in this work is to propose an appropriate MPC formulation to obtain a prescribed compliant behavior, while respecting a set of practical robotic constraints. The paper is organized as follows. In section II, we introduce the concepts of system dynamics and IC, later used in the paper. Section III describes the design method for the proposed MPC controller. Section IV provides experimental results featuring a collaborative robot. Finally, conclusions and perspectives are given in Section V.

II. CONTROL BACKGROUND A. Rigid Body Dynamics and Linearized Model

Consider the model of an n-joint serial robotic manipulator operating in an m-dimensional task space M (q)q + C( q, q) + G(q) = τ c -J(q) T f ext [START_REF] Goodrich | Human-robot interaction: A survey[END_REF] where q, q, q ∈ R n are joint accelerations, velocities and positions, respectively, τ c ∈ R n are the commanded joint torques, f ext ∈ R m is the end-effector wrench, J(q) ∈ R m×n is the robot Jacobian matrix and M (q) ∈ R n×n , C( q, q) ∈ R n and G(q) ∈ R n are the inertia matrix and the Coriolis and gravity terms, respectively. A classical solution [START_REF] Siciliano | Springer Handbook of Robotics[END_REF] to linearize equation ( 1) is to apply the control law τ c = M (q)ν + C( q, q) + G(q) + J(q) T f ext with ν the new control input, resulting in the double integrator model in joint space

q = ν (2)
The robot dynamics can then be expressed in task space by derivating the forward differential kinematic model, which leads to p = J(q)q + J(q, q) q (3)

where p ∈ R m , J(q, q) ∈ R m×n denote the robot endeffector pose and the time derivative of the Jacobian matrix. Using (2) and ( 3) one can express a new input u such that ν = J + u -J + J q, with J + ∈ R n×m the pseudo-inverse of J, resulting in a double integrator model in task space

p = u (4)

B. Impedance Control

In order to dynamically link the system positions, velocities and accelerations with the external forces, the input u of system (4) is computed in order to obtain

M v (p r -p) + D v ( ṗr -ṗ) + K v (p r -p) = f ext (5)
So, the system subject to the external forces f ext is characterized by its impedance with an apparent virtual mass M v ∈ R m×m , a desired damping D v ∈ R m×m and a desired stiffness K v ∈ R m×m , tracking a reference motion p r . Note that M v , D v and K v are symmetric positive definite matrices. The desired behavior (5) can be obtained by imposing

u = pr + M -1 v (D v ( ṗr -ṗ) + K v (p r -p) -f ext ) (6)
to the linearized system [START_REF] Siciliano | Springer Handbook of Robotics[END_REF]. IC can also be represented using a state space model. In the following, we will consider external forces with slow variations so that ḟext ≈ 0. This simplification corresponds to many practical cases but could seem restrictive in the case of unmodeled contact. However, if a better model of the external force is available, it can be included at this stage into the state-space model. This aspect of the problem goes far beyond the scope of this paper, justifying the previous assumption. With this assumption, the dynamic model of equation ( 4) can be written as

ẋ = A c x + B c u (7) with x =   ṗ p f ext   , A c =   0 0 0 I m 0 0 0 0 0   , B c =   I m 0 0  
The external force has been integrated into the state variable in order to better deal with this disturbance, inspired by the Internal Model Principle [START_REF] Francis | The internal model principle of control theory[END_REF]. Note however that the objective is not to cancel the disturbance, but rather to have an appropriate response with respect to it.

Let us now write the IC using a state space formalism. To do so, note that (6) can be written as

u = u r + K(x r -x) (8)
with

u r = pr ∈ R m (9)
x r = ṗT

r p T r 0 T ∈ R 3m (10) K = M -1 v D v K v I m ∈ R m×3m ( 11
)
where r is some reference vector that can be chosen according to the desired task. By substituting ( 8) into [START_REF] Kronander | Stability considerations for variable impedance control[END_REF], the state space model of the impedance controlled system can be written as

ẋ = (A c -B c K)x + B c (u r + Kx r ) (12)
IC can therefore be seen as a state-feedback tracking problem.

Since the controller is implemented numerically, it is necessary to find a discrete model of the system. The discretization of equation ( 7) with zero order hold yields

x k+1 = Ax k + Bu k (13)
where k represents the current step, x k is the discrete state at step k and A, B and C are the discrete state space matrices.

The sampling period T s is selected short enough to emulate the continuous controller [START_REF] Lasota | A survey of methods for safe humanrobot interaction[END_REF], which gives

u k = u r k -K(x k -x r k ) (14)
with x r k and u r k the discrete forms of x r and u r .

III. MODEL PREDICTIVE IMPEDANCE CONTROL

MPC consists in solving at every step a finite-horizon optimal control problem with constraints. In order to deal at the same time with unpredictable interactions as does the IC, and different types of limits as does the MPC, a novel controller, called Model Predictive Impedance Controller (MPIC) is proposed.

When no constraints are active, the MPIC should be equivalent to IC, i.e. to the static state feedback control [START_REF] Killpack | Model predictive control for fast reaching in clutter[END_REF]. Several methods have been proposed in the literature to make a MPC behave as a state feedback such as the LMIbased inverse optimality method [START_REF] Hartley | Designing output-feedback predictive controllers by reverse-engineering existing lti controllers[END_REF], [START_REF] Priess | Solutions to the Inverse LQR Problem With Application to Biological Systems Analysis[END_REF] or the controller matching based on QP matrices as in [START_REF] Cairano | Model predictive control tuning by controller matching[END_REF]. However, these methods either do not guarantee the exact matching between the MPC and the desired state feedback controller, or do not allow online computation. To obtain the controller of equation ( 14), the cost function J of the MPC problem has to be designed such that its optimum is zero when the desired behavior is obtained. To do so, a zero-value cost function over a control horizon H can be written, as in [START_REF] Hartley | Designing output-feedback predictive controllers by reverse-engineering existing lti controllers[END_REF].

J = H-1 k=0 (u k + K(x k -x r k ) -u r k ) T R(u k + K(x k -x r k ) -u r k ) (15)
which is optimal for the desired control law, for any R = R T > 0 ∈ R m×m . Furthermore, J can be written as

J = H-1 k=0 u k -u r k x k -x r k T R S T S Q u k -u r k x k -x r k ( 16
)
where

Q = K T RK ∈ R 3m×3m and S = K T R ∈ R 3m×m
are the weights of the particular terms satisfying the conditions Q = Q T ≥ 0 and Q -SR -1 S T ≥ 0, and H is the control horizon. In [START_REF] Kreindler | Optimality of linear control systems[END_REF] it is pointed out that the cross-terms between input and state in the quadratic cost function allow reproducing any state-feedback control strategy by making this equivalence as the primary control objective.

With the discrete state-space plant dynamics (13), a MPC problem can be written as the tracking of some references, with feedback x r k and u r k terms

min u k J = min u k H-1 k=0 u k -u r k x k -x r k T R S T S Q u k -u r k x k -x r k ( 17
)
s.t. x k+1 = Ax k + Bu k x 0 = measured current state
where x 0 is the initial state measured at each control cycle.

In the following, the input rate of change

∆u k = u k -u k-1
is used instead of the input u k , a choice often made for tracking problems. First, let us define

∆u =[∆u T 0 ... ∆u T H-1 ] T ∈ R mH u =[u T -1 u T 0 ... u T H-1 ] T ∈ R m(H+1) (18)
x T 0 ... x T H-1 ] T ∈ R 3mH ρ =[u r 0 T ... u r H-1 T x r 0 T ... x r H-1 T ] T ∈ R 4mH
the vectors containing all the input rates of change, inputs, states and references over the control horizon H. Notation u -1 represents the last applied input. The minimization problem [START_REF] Kazim | Combined predictive path following and admittance control[END_REF] can be transformed into a Quadratic Program (QP) form. To do so, it has first to be expressed in matrix form such that

min ∆u J = min ∆u     ∆u x u ρ     T Φ T R S T S Q Φ     ∆u x u ρ     (19)
s.t. x = Ax 0 + Cu u = T ∆u + Iu -1 (20)
with u, ∆u, x and ρ defined by [START_REF] Francis | The internal model principle of control theory[END_REF] and

Q =diag(Q, ... , Q) ∈ R 3mH×3mH R =diag(R, ... , R) ∈ R mH×mH S =diag(S, ... , S) ∈ R 3mH×mH I =[I m ... I m ] T ∈ R m(H+1)×m
and block matrices H+1) is the H-steps state reachability matrix and A ∈ R 3mH×3m the H-steps free evolution matrix [START_REF] Cairano | Model predictive control tuning by controller matching[END_REF], such that 19) is a constrained QP that could be solved numerically. However, previously, a simplification of the problem can be performed by including the constraints (20) into the cost function, so that

T =         0 0 . . . 0 I m 0 . . . 0 . . . . . . . . . . . . I m . . . . . . 0 I m . . . . . . I m         ∈ R m(H+1)×mH Φ = I mH 0 I mH 0 -I mH 0 0 I 3mH 0 0 0 -I 3mH ∈ R 4mH×m(5H+1) C ∈ R 3mH×m(
C =        0 0 0 . . . 0 0 B 0 . . . 0 0 AB B . . . 0 . . . . . . . . . . . . . . . 0 A H-2 B A H-3 B . . . B        A =        I 3m A A 2 . . . A H-1        Problem (
J =     ∆u x 0 u -1 ρ     T Ψ T Φ T R S T S Q ΦΨ     ∆u x 0 u -1 ρ     (21)
with

Ψ =     I mH 0 0 0 CT A CI 0 T 0 I 0 0 0 0 I 4mH     ∈ R m(5H+1)×m(5H+4)
By defining matrices

H ∈ R mH×mH , F ∈ R mH×m(5H+4) and L ∈ R m(5H+4)×m(5H+4) such that H F F T L = Ψ T Φ T R S T S Q ΦΨ ( 22
)
the optimization problem [START_REF] Hartley | Designing output-feedback predictive controllers by reverse-engineering existing lti controllers[END_REF] becomes

min ∆u J = min ∆u ∆u T H∆u + 2[x T 0 u T -1 ρ T ]F T ∆u (23)
where the term corresponding to L has been removed as it is constant. Because of the hypotheses on the weight matrices, H > 0 and the equation ( 23) is a convex QP. Note here, that as the problem is convex the solution corresponding to the zero cost function is the unique global solution. This means, that when unconstrained, the optimal input given by the MPIC controller will behave as the desired state-feedback. Without constraints, [START_REF] Chawda | Toward controlling a kuka lbr iiwa for interactive tracking[END_REF] has the following analytic solution

∆u = -H -1 F x T 0 u T -1 ρ T T (24)
In MPC, only the first input is applied at each iteration, i.e. u 0 = u -1 -I m 0 ... 0 ∆u. Note that this applies only when there are no other constraints than the system model. However, in most practical cases, additional constraints have to be considered. In this case, no analytical solution exists and a numerical solver has to be used. In order to do so, the constraints need to be written in the standard form

Γ∆u ≤ Π + Ω x T 0 u T -1 T ( 25
)
where Γ ∈ R c×mH , Ω ∈ R c×4m and Π ∈ R c are formulated based on c problem constraints.

IV. EVALUATION A. Experimental setup

The performance of the proposed controller was evaluated on an experimental setup featuring a collaborative robot following a reference trajectory. We chose to perform the experiments on a simple exemplary environment in order to ensure good repeatability of the task, and so to provide comparable results. For this reason, no interaction with the user occurs during the experiments presented in the paper. However, human-robot interactions are shown in the accompanying video. The reference was a smooth polynomial trajectory on a planar surface, placed slightly under the table in order to ensure contact of the end-effector with the (x,y)plane. In addition, a rigid obstacle was placed on the path generating unmodeled disturbances. The experimental setup is shown in Fig. 1. Both IC and MPIC were compared, throughout four different experiments. At first, we show the equivalence of IC and MPIC when no constrains are applied. Then, position, velocity and input constraints are successively imposed to the MPIC controller for validation of its performances, and to illustrate the advantages it offers.

B. Implementation and hardware specifications

The experiments were carried out with a KUKA iiwa 14 collaborative robot. An ATI-Mini40 force-torque sensor was mounted on the robot to measure the end-effector wrench. The proposed controller addresses the robot at joint torque level, which is not a classical approach for controlling this type of manipulators. As the internal structure of the low level controllers of the KUKA iiwa are not well known, with some identification attempts such as in [START_REF] Chawda | Toward controlling a kuka lbr iiwa for interactive tracking[END_REF], [START_REF]Toward torque control of a kuka lbr iiwa for physical humanrobot interaction[END_REF], we based our understanding of the low level control architecture on [START_REF] Albu-Schffer | A unified passivity-based control framework for position, torque and impedance control of flexible joint robots[END_REF]. The robot software allows the tuning of joint impedance parameters, which were set to zero and the torque commands generated by the MPIC were fed to the system using KUKA's FRI (Fast Robot Interface) protocol. This protocol was also used to measure the joint positions, whereas the joint velocities were computed using a filtered derivative with a cut-off frequency at 50Hz. End-effector positions and velocities were computed using the forward kinematic and forward differential kinematic models, respectively. The control and the FRI communication were sampled at 500Hz. For solving the QP problem [START_REF] Chawda | Toward controlling a kuka lbr iiwa for interactive tracking[END_REF] with constraints, the C++ solver qpOASES [START_REF] Ferreau | qpoases: a parametric active-set algorithm for quadratic programming[END_REF] was used. The impedance parameters were chosen such that M v = diag(5, 5, 5, 2, 2, 2), K v = diag(300, 300, 50, 500, 500, 500), and D v = diag(77.5, 77.5, 31.6, 63.2, 63.2, 63.2). K v was chosen in order to ensure: 1) a good contact in the direction normal to the (x,y)-plane, 2) a relatively low stiffness in the x and y directions, which represents a trade-off in compliance and tracking performance, and 3) a high stiffness in orientation in order to have a nearly constant wrist pose. D v was selected in order to have a damping ratio of 1 for the impedance model. The prediction horizon for the MPIC was set to H = 5, as currently the solver struggles to solve the constrained QP problem for a longer horizon in the given time. As the input given by the MPIC when operating close to a constraint can be noisy, a fixed limit on the input rate of change -0.1 ≤ ∆u ≤ 0.1 was set.

C. Experimental results

1) Comparing IC and MPIC without constraints: The IC law ( 6) was implemented on the system and tested in the same conditions as the MPIC. Fig. 2 allows comparing the two controllers. At position level (Fig. 2a) one can observe the almost perfect matching of the two trajectories. The position error between IC and MPIC has been evaluated to ēx = 0.24mm and ēy = 0.12mm with a standard deviation σ x = 1.3mm and σ y = 0.7mm. This match can also be noticed at velocity (Fig. 2b) and acceleration (Fig. 2c) levels where only a small difference on the x-axis can be spotted.

We believe that this difference is caused by some friction on the surface. One can also notice that both controllers do not perfectly match the reference trajectory. The main reason lies in the fact that the impedance model has a rather low stiffness, which impacts reference tracking performance. However, at some points, as it is the case in the left upper corner of the trajectory in Fig. 2a, the performance of both controllers seems to be particularly reduced, even though they behave identically. We believe that this can be caused by the internal controller of the KUKA robot that has some issues tracking torque commands as pointed out in [START_REF] Chawda | Toward controlling a kuka lbr iiwa for interactive tracking[END_REF], [START_REF]Toward torque control of a kuka lbr iiwa for physical humanrobot interaction[END_REF].

We were interested in verifying the impact of these small defects due to the (not open) low level controllers of KUKA iiwa. To do so, we identified the real impedance parameters estimated from experimental data, and compared them to the ones chosen in the impedance model. A linear least square minimization was performed on the experimental data of IC to retrieve these parameters. The resulting estimations are given in Table I. It turns out that the system impedance during the experiments fits rather well the imposed impedance model, with a 8% max error on M v , a 6.6% max error on K v , and a 6.1% max error on D v , for the x and y components. This result shows that the uncertainties on the internal controller do not heavily impact the performance and do not interfere with the robot-environment interaction. Typically, a human operator will not be able to discriminate such small differences in the inertia, damping and stiffness properties [START_REF] Jones | Kinesthetic sensing[END_REF].

2) Position constrained MPIC: In this experiment, position constraints were specified to the MPIC. Fig. 3 shows the resulting end-effector trajectories as well as the position constraint, which was set such that p x ≥ -0.006m. One can observe that the two controllers have the same behavior when working in the allowed task space. When approaching the limit, MPIC stops behaving as IC to meet the constraint. End-effector trajectories. Case 2: IC and MPIC with position constraints (x ≥ -0.006m) and with an obstacle on the path.

3) Velocity constrained MPIC: We show next, that it is also possible to apply velocity constraints with the MPIC controller. Fig. 4 compares the behavior of MPIC and IC at position and velocity level. In Fig. 4b, one can see that given velocity constraints v x ≤ 0.02m/s, v y ≤ 0.02m/s, the MPIC is able to limit the system velocity below the desired threshold. As in the previous paragraph, this also results in the modification of the desired trajectory, as shown Fig. 4a. 4) Input constrained MPIC: Finally, the MPIC was tested in the case of input constraints resulting in the limitation of the commanded acceleration of the system. Fig. 5a allows comparing the accelerations obtained via IC and MPIC for a constraint a y ≥ -0.4m/s 2 . As observed previously, both controllers have an equivalent behavior when inside the allowed acceleration range and MPIC limits acceleration when exceeding the limit. It is worth noticing here that limiting the acceleration will affect the contact force between the robot and the obstacle, as shown in Fig. 5b. This feature is of particular interest as it acts as a saturation that can limit the input to some maximally allowed threshold. In addition, it also acts directly on the system maximum acceleration and the resulting force applied by the robot, which is an important feature for interaction tasks. Experimental comparison of acceleration and measured force applied by the robotic manipulator end-effector. Case 4: IC and MPIC control with acceleration constraints (ay ≥ -0.4m/s 2 ) and with an obstacle on the path.

V. CONCLUSION

In this paper, we have proposed a method to synthesize an impedance controller (IC) as a model predictive controller (MPC), based on the state-space representation of IC. This new controller is called Model Predictive Impedance Controller (MPIC). We have shown that MPIC combines the features of IC with the advantages of MPC, allowing to impose constraints on different system parameters, such as state or input. The controller has been tested with an experimental setup featuring a collaborative robot, and compared to IC. MPIC behaves exactly as IC when no constraints are active and it respects position, velocity and acceleration constraints when such constraints are imposed.

These results open the way to many possible applications, in particular to ensure safety in critical applications of col-laborative robotics. Future work will tend towards including a human operator into the interaction dynamics to take full advantage of the prediction capability of MPIC to anticipate interactions. Another enhancement will be to improve the current setup, so that MPIC might be capable of adapting the impedance parameters online during the execution of the task. In addition, a theoretical stability analysis will be investigated in order to specify possible instability issues of MPIC, as currently only local stability can be guaranteed. One way to do so would be to impose some additional constraints that limit the system global energy, or to mathematically prove that when operating close to boundaries the system remains stable.

Fig. 1 .

 1 Fig. 1. Experimental setup with reference trajectory and obstacle.

  (a) End-effector trajectories under IC and MPIC. (b) End-effector velocities under IC and MPIC. (c) End-effector accelerations under IC and MPIC.

Fig. 2 .

 2 Fig. 2.Experimental comparison of the trajectories, velocities and accelerations of the robotic manipulator end-effector, following a reference trajectory. Case 1: IC and MPIC without active constraints and with an obstacle on the path.

Fig. 3 .

 3 Fig. 3.End-effector trajectories. Case 2: IC and MPIC with position constraints (x ≥ -0.006m) and with an obstacle on the path.

  (a) End-effector trajectories under IC and MPIC. (b) End-effector velocities under IC and MPIC.

Fig. 4 .

 4 Fig. 4. Experimental comparison of the trajectories and velocities of the robotic manipulator end-effector following a reference trajectory. Case 3: IC and MPIC with active velocity constraints (vx ≤ 0.02m/s, vy ≤ 0.02m/s) and with an obstacle on the path.

  (a) End-effector accelerations under IC and MPIC. (b) Force measured at the end-effector under IC and MPIC.

Fig. 5 .

 5 Fig. 5.Experimental comparison of acceleration and measured force applied by the robotic manipulator end-effector. Case 4: IC and MPIC control with acceleration constraints (ay ≥ -0.4m/s 2 ) and with an obstacle on the path.

TABLE I IMPEDANCE

 I PARAMETERS RETRIEVED FROM DATA.

		Mv	Kv	Dv
	x	5.40 290.04 76.72
	y	5.39 280.27 72.80

ACKNOWLEDGMENT This work was supported by the Investissements d'Avenir program of the French Government, more specifically by Robotex ANR-10-EQPX-44, Labex CAMI ANR-11-LABX-0004 and the IDEX of the University of Strasbourg.