We investigate the numerical solution of stable Sylvester equations via iterative schemes proposed for computing the sign function of a matrix. In particular, we discuss how the rational iterations for the matrix sign function can efficiently be adapted to the special structure implied by the Sylvester equation. For Sylvester equations with factored constant term as those arising in model reduction or image restoration, we derive an algorithm that computes the solution in factored form directly. We also suggest convergence criteria for the resulting iterations and compare the accuracy and performance of the resulting methods with existing Sylvester solvers. The algorithms proposed here are easy to parallelize. We report on the parallelization of those algorithms and demonstrate their high efficiency and scalability using experimental results obtained on a cluster of Intel Pentium Xeon processors.

We investigate the numerical solution of stable Sylvester equations via iterative schemes proposed for computing the sign function of a matrix. In particular, we discuss how the rational iterations for the matrix sign function can efficiently be adapted to the special structure implied by the Sylvester equation. For Sylvester equations with factored constant term as those arising in model reduction or image restoration, we derive an algorithm that computes the solution in factored form directly. We also suggest convergence criteria for the resulting iterations and compare the accuracy and performance of the resulting methods with existing Sylvester solvers. The algorithms proposed here are easy to parallelize. We report on the parallelization of those algorithms and demonstrate their high efficiency and scalability using experimental results obtained on a cluster of Intel Pentium Xeon processors.

Given a full-range simply-invariant shift-invariant subspace <i>M</i> of the vector-valued <i>L²</i> space on the unit circle, the classical Beurling-Lax-Halmos (BLH) theorem obtains a unitary operator-valued function <i>W</i> so that <i>M</i> may be represented as the image of of the Hardy space <i>H²</i> on the disc under multiplication by <i>W</i>. The work of Ball-Helton later extended this result to find a single function representing a so-called dual shift-invariant pair of subspaces <i>(M,M^Ã)</i> which together form a direct-sum decomposition of <i>L²</i>. In the case where the pair <i>(M,M^Ã)</i> are finite-dimensional perturbations of the Hardy space <i>H²</i> and its orthogonal complement, Ball-Gohberg-Rodman obtained a transfer function realization for the representing function <i>W</i>; this realization was parameterized in terms of zero-pole data computed from the pair <i>(M,M^Ã)</i>. Later work by Ball-Raney extended this analysis to the case of nonrational functions <i>W</i> where the zero-pole data is taken in an infinite-dimensional operator theoretic sense. The current work obtains analogues of these various results for arbitrary dual shift-invariant pairs <i>(M,M^Ã)</i> of the <i>L²</i> spaces on the real line; here, shift-invariance refers to invariance under the translation group. These new results rely on recent advances in the understanding of continuous-time infinite-dimensional input-state-output linear systems which have been codified in the book by Staffans. / Ph. D.

L’imagerie hyperspectrale (HS) consiste à acquérir une même scène dans plusieurs centaines de bandes spectrales contiguës (dimensions d'un cube de données), ce qui a conduit à trois types d'applications pertinentes, telles que la détection de cibles, la classification et le démélange spectral. Cependant, tandis que les capteurs hyperspectraux fournissent une information spectrale abondante, leur résolution spatiale est généralement plus limitée. Ainsi, la fusion d’une image HS avec d'autres images à haute résolution de la même scène, telles que les images multispectrales (MS) ou panchromatiques (PAN) est un problème intéressant. Le problème de fusionner une image HS de haute résolution spectrale mais de résolution spatiale limitée avec une image auxiliaire de haute résolution spatiale mais de résolution spectrale plus limitée (parfois qualifiée de fusion multi-résolution) a été exploré depuis de nombreuses années. D'un point de vue applicatif, ce problème est également important et est motivé par ceratins projets, comme par exemple le project Japonais HISIU, qui vise à fusionner des images MS et HS recalées acquises pour la même scène avec les mêmes conditions. Les techniques de fusion bayésienne permettent une interprétation intuitive du processus de fusion via la définition de la loi a posteriori de l’image à estimer (qui est de hautes résolutions spatiale et spectrale). Puisque le problème de fusion est généralement mal posé, l’inférence bayésienne offre un moyen pratique pour régulariser le problème en définissant une loi a priori adaptée à la scène d'intérêt. Les différents chapitres de cette thèse sont résumés ci-dessous. Le introduction présente le modèle général de fusion et les hypothèses statistiques utilisées pour les images multi-bandes observées, c’est-à-dire les images HS, MS ou PAN. Les images observées sont des versions dégradées de l'image de référence (à hautes résolutions spatiale et spectrale) qui résultent par exemple d’un flou spatial et spectral et/ou d’un sous-échantillonnage liés aux caractéristiques des capteurs. Les propriétés statistiques des mesures sont alors obtenues directement à partir d’un modèle linéaire traduisant ces dégradations et des propriétés statistiques du bruit. Le chapitre 1 s’intéresse à une technique de fusion bayésienne pour les images multi-bandes de télédétection, à savoir pour les images HS, MS et PAN. Tout d'abord, le problème de fusion est formulé dans un cadre d'estimation bayésienne. Une loi a priori Gaussienne exploitant la géométrie du problème est définie et un algorithme d’estimation Bayésienne permettant d’estimer l’image de référence est étudié. Pour obtenir des estimateurs Bayésiens liés à la distribution postérieure résultant, deux algorithmes basés sur échantillonnage de Monte Carlo et l'optimisation stratégie ont été développés. Le chapitre 2 propose une approche variationnelle pour la fusion d’images HS et MS. Le problème de fusion est formulé comme un problème inverse dont la solution est l'image d’intérêt qui est supposée vivre dans un espace de dimension résuite. Un terme de régularisation imposant des contraintes de parcimonie est défini avec soin. Ce terme traduit le fait que les patches de l'image cible sont bien représentés par une combinaison linéaire d’atomes appartenant à un dictionnaire approprié. Les atomes de ce dictionnaire et le support des coefficients des décompositions des patches sur ces atomes sont appris à l’aide de l’image de haute résolution spatiale. Puis, conditionnellement à ces dictionnaires et à ces supports, le problème de fusion est résolu à l’aide d’un algorithme d’optimisation alternée (utilisant l’algorithme ADMM) qui estime de manière itérative l’image d’intérêt et les coefficients de décomposition. / Hyperspectral (HS) imaging, which consists of acquiring a same scene in several hundreds of contiguous spectral bands (a three dimensional data cube), has opened a new range of relevant applications, such as target detection [MS02], classification [C.-03] and spectral unmixing [BDPD+12]. However, while HS sensors provide abundant spectral information, their spatial resolution is generally more limited. Thus, fusing the HS image with other highly resolved images of the same scene, such as multispectral (MS) or panchromatic (PAN) images is an interesting problem. The problem of fusing a high spectral and low spatial resolution image with an auxiliary image of higher spatial but lower spectral resolution, also known as multi-resolution image fusion, has been explored for many years [AMV+11]. From an application point of view, this problem is also important as motivated by recent national programs, e.g., the Japanese next-generation space-borne hyperspectral image suite (HISUI), which fuses co-registered MS and HS images acquired over the same scene under the same conditions [YI13]. Bayesian fusion allows for an intuitive interpretation of the fusion process via the posterior distribution. Since the fusion problem is usually ill-posed, the Bayesian methodology offers a convenient way to regularize the problem by defining appropriate prior distribution for the scene of interest. The aim of this thesis is to study new multi-band image fusion algorithms to enhance the resolution of hyperspectral image. In the first chapter, a hierarchical Bayesian framework is proposed for multi-band image fusion by incorporating forward model, statistical assumptions and Gaussian prior for the target image to be restored. To derive Bayesian estimators associated with the resulting posterior distribution, two algorithms based on Monte Carlo sampling and optimization strategy have been developed. In the second chapter, a sparse regularization using dictionaries learned from the observed images is introduced as an alternative of the naive Gaussian prior proposed in Chapter 1. instead of Gaussian prior is introduced to regularize the ill-posed problem. Identifying the supports jointly with the dictionaries circumvented the difficulty inherent to sparse coding. To minimize the target function, an alternate optimization algorithm has been designed, which accelerates the fusion process magnificently comparing with the simulation-based method. In the third chapter, by exploiting intrinsic properties of the blurring and downsampling matrices, a much more efficient fusion method is proposed thanks to a closed-form solution for the Sylvester matrix equation associated with maximizing the likelihood. The proposed solution can be embedded into an alternating direction method of multipliers or a block coordinate descent method to incorporate different priors or hyper-priors for the fusion problem, allowing for Bayesian estimators. In the last chapter, a joint multi-band image fusion and unmixing scheme is proposed by combining the well admitted linear spectral mixture model and the forward model. The joint fusion and unmixing problem is solved in an alternating optimization framework, mainly consisting of solving a Sylvester equation and projecting onto a simplex resulting from the non-negativity and sum-to-one constraints. The simulation results conducted on synthetic and semi-synthetic images illustrate the advantages of the developed Bayesian estimators, both qualitatively and quantitatively.

Nous nous intéressons dans cette thèse, à l’étude des méthodes itératives pour la résolutiond’équations matricielles de grande taille : Lyapunov, Sylvester, Riccati et Riccatinon symétrique.L’objectif est de chercher des méthodes itératives plus efficaces et plus rapides pour résoudreles équations matricielles de grande taille. Nous proposons des méthodes itérativesde type projection sur des sous espaces de Krylov par blocs Km(A, V ) = Image{V,AV, . . . ,Am−1V }, ou des sous espaces de Krylov étendus par blocs Kem(A, V ) = Image{V,A−1V,AV,A−2V,A2V, · · · ,Am−1V,A−m+1V } . Ces méthodes sont généralement plus efficaces et rapides pour les problèmes de grande dimension. Nous avons traité d'abord la résolution numérique des équations matricielles linéaires : Lyapunov, Sylvester, Stein. Nous avons proposé une nouvelle méthode itérative basée sur la minimisation de résidu MR et la projection sur des sous espaces de Krylov étendus par blocs Kem(A, V ). L'algorithme d'Arnoldi étendu par blocs permet de donner un problème de minimisation projeté de petite taille. Le problème de minimisation de taille réduit est résolu par différentes méthodes directes ou itératives. Nous avons présenté ainsi la méthode de minimisation de résidu basée sur l'approche global à la place de l'approche bloc. Nous projetons sur des sous espaces de Krylov étendus Global Kem(A, V ) = sev{V,A−1V,AV,A−2V,A2V, · · · ,Am−1V,A−m+1V }. Nous nous sommes intéressés en deuxième lieu à des équations matricielles non linéaires, et tout particulièrement l'équation matricielle de Riccati dans le cas continu et dans le cas non symétrique appliquée dans les problèmes de transport. Nous avons utilisé la méthode de Newtown et l'algorithme MINRES pour résoudre le problème de minimisation projeté. Enfin, nous avons proposé deux nouvelles méthodes itératives pour résoudre les équations de Riccati non symétriques de grande taille : la première basée sur l'algorithme d'Arnoldi étendu par bloc et la condition d'orthogonalité de Galerkin, la deuxième est de type Newton-Krylov, basée sur la méthode de Newton et la résolution d'une équation de Sylvester de grande taille par une méthode de type Krylov par blocs. Pour toutes ces méthodes, les approximations sont données sous la forme factorisée, ce qui nous permet d'économiser la place mémoire en programmation. Nous avons donné des exemples numériques qui montrent bien l'efficacité des méthodes proposées dans le cas de grandes tailles. / In this thesis, we focus in the studying of some iterative methods for solving large matrix equations such as Lyapunov, Sylvester, Riccati and nonsymmetric algebraic Riccati equation. We look for the most efficient and faster iterative methods for solving large matrix equations. We propose iterative methods such as projection on block Krylov subspaces Km(A, V ) = Range{V,AV, . . . ,Am−1V }, or block extended Krylov subspaces Kem(A, V ) = Range{V,A−1V,AV,A−2V,A2V, · · · ,Am−1V,A−m+1V }. These methods are generally most efficient and faster for large problems. We first treat the numerical solution of the following linear matrix equations : Lyapunov, Sylvester and Stein matrix equations. We have proposed a new iterative method based on Minimal Residual MR and projection on block extended Krylov subspaces Kem(A, V ). The extended block Arnoldi algorithm gives a projected minimization problem of small size. The reduced size of the minimization problem is solved by direct or iterative methods. We also introduced the Minimal Residual method based on the global approach instead of the block approach. We projected on the global extended Krylov subspace Kem(A, V ) = Span{V,A−1V,AV,A−2V,A2V, · · · ,Am−1V,A−m+1V }. Secondly, we focus on nonlinear matrix equations, especially the matrix Riccati equation in the continuous case and the nonsymmetric case applied in transportation problems. We used the Newton method and MINRES algorithm to solve the projected minimization problem. Finally, we proposed two new iterative methods for solving large nonsymmetric Riccati equation : the first based on the algorithm of extended block Arnoldi and Galerkin condition, the second type is Newton-Krylov, based on Newton’s method and the resolution of the large matrix Sylvester equation by using block Krylov method. For all these methods, approximations are given in low rank form, wich allow us to save memory space. We have given numerical examples that show the effectiveness of the methods proposed in the case of large sizes.