50 Matemática Universitária nº46

Artigo

Estimador MUSIC: resolvendoum problema de probabilidadecom álgebra linearPhilip Thompson1

IMPA

U m problema geral na teoria da informação (queenvolve, por exemplo, a Internet) consiste em

transmitir um sinal, que pode ser uma sequência finita(Fk)m

k=1, enviado de um emissor para um receptor atra-vés de um canal de comunicação. Um canal de comuni-cação pode significar fisicamente diferentes meios: umcabo coaxial de uma televisão, um fio de uma linha detransmissão da rede elétrica, a atmosfera por onde sepropaga um sinal de satélite até a antena de uma casaou celular ou, até mesmo, o ambiente de uma sala noqual uma rede wireless é instalada e um pulso eletro-magnético se propaga de um roteador até um notebook.Num canal não trivial, o sinal (Fk)m

k=1 é alterado indese-javelmente e o receptor recebe um sinal distinto, isto é,uma sequência finita (Xk)m

k=1. Mais precisamente, poralteração do sinal (Fk)m

k=1 queremos dizer que o sinalXk, que chega ao receptor, é o sinal Fk emitido maisum sinal indesejável Wk, que surge inevitavelmente du-rante a transmissão pelo canal. Isto é, Xk = Fk + Wk,para k = 1, . . . , m. É comum designar o sinal (Wk)m

k=1por ruído. Dependendo do modelo físico do canal, oruído surge devido a vários fatores. Por exemplo, umafonte comum de ruído num cabo coaxial é o ruído tér-mico decorrente da agitação térmica dos elétrons.

O ruído é um sinal indesejado que se manisfesta adi-tivamente e aleatoriamente. O principal objetivo da teo-ria da informação é recuperar, da maneira mais fiel pos-sível, o sinal (Fk)m

k=1 enviado e corrompido pelo ruído,conhecendo-se o sinal (Xk)m

k=1 medido no receptor. Ateoria da probabilidade se constitui numa ferramentaútil neste tipo de problema pois define rigorosamente

1 O autor e mestre pelo IMPA.

a característica aleatória do ruído. Mais precisamente,o ruído é matematicamente modelado como uma variá-vel aleatória, conceito bem definido em probabilidade, oqual veremos adiante.

Em 1979, R. O. Schmidt ([7]) elaborou um algoritmomuito eficiente capaz de solucionar o problema descritoagora. Considere p senoides complexas discretas

γ1eiθ1k, . . . , γpeiθpk, k ∈ N

com amplitudes γ1, . . . , γp e frequências θ1, . . . , θp.Considere também uma sequência (Wk)k∈N, que éum ruído gaussiano branco (cuja definição será vistamais adiante). A pergunta consiste em determinar asfrequências θ1, . . . , θp dessas senoides tendo informaçãoestatística da soma das senoides mais o ruído após munidades de tempo, isto é, da sequência finita (Xk)m

k=1dada por

Xk =

p

∑j=1

γjeiθjk + Wk : k = 1, . . . , m

.

Para facilitar a referência posterior, chamaremos esseproblema de o problema das senoides. Ele é um casoparticular do problema de transmissão por um canalde comunicação descrito anteriormente, em que o si-nal transmitido é dado por Fk = ∑

pj=1 γje

iθjk, isto é, umsinal senoidal (soma de senoides complexas). Convémcitar que, sob a ação aleatória do ruído, as amplitudesγ1, . . . , γp das senoides são também modeladas comovariáveis aleatórias. Note que, sob influência do ruído,a sequência medida X := (Xk)m

k=1 também é uma va-riável aleatória e veremos que as informações estatís-ticas suficientes referidas acima serão sua média e suacovariância. Apesar de definidas abstratamente (o queserá feito adiante), essas quantidades podem ser apro-ximadamente estimadas com precisão desde que o nú-mero de medições seja suficientemente grande. Esta éa filosofia da estatística matemática, que torna a teoria

1

51Matemática Universitária nº46

Artigo

probabilística aplicável a experimentos concretos. Porexemplo: se, durante o experimento, são medidos K va-lores X(ωi) : i = 1, . . . , K da variável aleatória X refe-rentes aos resultados ωi; i = 1, . . . , K do experimento,então a estimativa mais comum para a média de X é

E(X) :=1K

K

∑i=1

X(ωi) ,

isto é, a média aritmética entre os valores medidos. Aestimativa da covariância segue um formato análogo.

A formulação rigorosa desse problema no contextoda teoria da probabilidade (o que tornará preciso o sig-nificado de aleatoriedade descrita acima) será dada adi-ante (na Seção 2). Historicamente, a motivação para oproblema surgiu num contexto de aplicação em proces-samento de sinais (veja [7]). Neste sentido, um dos ob-jetivos deste artigo é apresentar conceitos de probabili-dade ao estudante de matemática, visando formalizarum problema motivado pelas aplicações. Mais do queisso, o objetivo deste artigo é apresentar um exemplode um problema de probabilidade cuja solução requeressencialmente conceitos e resultados da álgebra linear(tais como operadores autoadjuntos e o Teorema Espec-tral). Neste sentido, o artigo é acessível a estudantes degraduação e constrói uma ponte entre uma disciplinabásica e uma área mais avançada. O algoritmo é deno-minado MUSIC (MUltiple SIgnal Classification). Esteartigo apresenta o problema e sua solução.

Para motivar a leitura, é importante citar que umgrande número de aplicações na área de comunicaçõespode ser modelado pelo problema das senoides citadoacima. Dentre elas, podemos citar: processamento desinais de áudio, voz e imagem; estimativa da direçãode chegada em radares e sonares; e ressonância nuclearmagnética. Talvez o contexto mais natural de aplicaçãocujo modelo recai no problema das senoides surja nosistema de modulação em frequência utilizado na trans-missão e recepção de sinais enviados por radar (veja [6]e [1]). A utilidade de um radar consiste na possibilidadede detectar a presença, distância e velocidade de um ob-jeto ao alcance de medição. Sabe-se, em comunicações,que devido às limitações físicas do canal de comuni-cação (tais como atenuação do sinal), torna-se inviáveltransmitir a informação diretamente. Uma estratégiacomum para superar tal problema é a chamada modu-

lação. Tal técnica consiste em codificar o sinal transmi-tido num outro sinal mais adequado para a transmis-são, por exemplo um sinal senoidal. Neste caso, o si-nal a ser transmitido pode ser codificado, por exemplo,na frequência do sinal senoidal (este tipo de modula-ção é chamada de modulação em frequência). Grosso modo,um pulso eletromagnético senoidal é emitido pelo ra-dar, que se propaga até atingir o alvo. Após atingi-lo,o pulso é refletido de volta ao radar e medido. Pode-seusar técnicas apropriadas de modulação de sinais pro-jetadas no circuito do radar, de tal maneira que o pulsosenoidal recebido por este radar possua frequência pro-porcional à distância do alvo ao radar. A constante deproporcionalidade só depende da velocidade de propa-gação do pulso (isto é, a velocidade da luz) e o tempoexpirado entre a transmissão e a recepção do pulso. Ouseja, descobrindo a frequência do sinal medido, desco-brimos a distância do alvo. A rigor, isto funciona bem seo alvo está em repouso. No caso em que o alvo está emmovimento, outras técnicas são necessárias para incor-porar o chamado Efeito Doppler, no qual a frequênciadepende da velocidade do alvo relativamente ao radar.

Apesar de não detalharmos aqui, convém citar o con-texto original de aplicação em que o algoritmo MUSICfoi proposto por Schmidt, em [7]. Usando-se o mesmoprincípio de transmissão e recepção de sinais via ra-dar, MUSIC pode ser utilizado para estimar a orien-tação espacial de um objeto no espaço bi ou tridimen-sional. Mais especificadamente, MUSIC permite esti-mar o ângulo-de-chegada do pulso eletromagnético re-fletido pelo alvo até o radar. Neste caso, o modelo nãoé dado pelo problema das senoides, mas como menci-onaremos, MUSIC pode ser usado num modelo linearmais geral.

Outra aplicação de MUSIC surge em processamentode imagens obtidas por radar. Através de um esquemacomplicado de processamento, imagens bidimensio-nais e tridimensionais de uma vasta região podem serproduzidas utilizando-se o mesmo princípio de trans-missão e recepção de pulsos eletromagnéticos, via ra-dares instalados em aviões e satélites, que percorrema região emitindo pulsos e medindo os refletidos pelosolo repetidas vezes. Neste contexto, o problema das se-noides também surge e MUSIC pode ser utilizado (veja

2

52 Matemática Universitária nº46

Artigo

[2]).

Outro contexto comum em que o parâmetro a ser es-timado é codificado na frequência de uma senoide é ode processamento de áudio. Como se sabe da análiseharmônica, sob certas condições um sinal contínuo éessencialmente uma soma de várias senoides compo-nentes. O conjunto destas frequências é denominadoespectro do sinal. Em processamento de áudio e voz,muitas vezes deseja-se detectar uma frequência princi-pal de um sinal imerso em ruído. Vale a pena citar quea Transformada de Fourier (FT) é uma ferramenta pa-drão nesse tipo de problema em processamento de si-nais. MUSIC revelou-se como um método alternativoa FT por obter uma resolução significativamente maior.Sua desvantagem em relação à FT é a necessidade deuma amostra estatística do sinal, o que torna seu usoinviável em contextos em que não está disponível umaamostra suficiente. Entretanto, em casos nos quais exis-tem poucas amostras disponíveis, mas certas hipóte-ses de ergodicidade do sinal são satisfeitas, o algoritmoMUSIC pode ser usado conjuntamente com o métodospatial smoothing, em que partes da amostra são trata-das como diferentes amostras do sinal, aumentando onúmero de amostras disponíveis ([2]).

Enfim, o problema das senoides é um exemplo de pro-blema inverso estatístico, que é como se formulam ma-tematicamente variados problemas de engenharia, eco-nomia, física e outras áreas da ciência.

Nesta primeira seção, apresentamos de forma sucintadefinições fundamentais em probabilidade que nos for-necem a linguagem básica no estudo de sistemas pro-babilísticos. Não é assumido que o leitor esteja familia-rizado com conceitos de teoria da medida. Na verdade,não usaremos tecnicamente, nas provas a seguir, qual-quer resultado dessa teoria. Entretanto, precisamos sa-ber o que é uma variável aleatória para entendermos ma-tematicamente o que vem a ser “ruído”. Além disso,precisamos saber o que é valor esperado pois ele permi-tirá definir a média e a variância de uma variável alea-tória que, grosso modo, caracterizam as variáveis aleató-rias. A apresentação a seguir será devidamente resu-

mida de forma a atender esses propósitos. O leitor in-teressado numa apresentação mais completa pode con-sultar os primeiros capítulos de [3]. Nesta breve expo-sição, apresentaremos o que vem a ser espaço de probabi-lidade, variável aleatória, valor esperado, variância e covari-ância.

1.1 Espaco de probabilidade

Temos disponível um conjunto Ω e uma família M desubconjuntos de Ω. Em probabilidade, interpreta-se Ωcomo o conjunto de todos os resultados possíveis de umexperimento, e cada elemento de M é denominado umevento do experimento. O conjunto Ω é denominado es-paço amostral. Por exemplo, se o experimento em ques-tão é o arremesso de um dado, o espaço amostral éΩ = 1, 2, 3, 4, 5, 6 e M é o conjunto das partes de Ω.O evento “tirar faces 1 ou 2” é o subconjunto 1, 2. Naverdade, a família M não pode ser qualquer e exige-seque ela permita operar subconjuntos de Ω pelas opera-ções de união, intersecção e complementaridade. Maisprecisamente, exige-se que:

• Ω ∈ M;

• se A ∈ M então AC ∈ M;

• se (Ai)∞i=1 é uma família enumerável de elementos

de M então∞

i=1 Ai ∈ M.

Uma família M com tais propriedades é dita uma σ-álgebra de Ω. Cada subconjunto de Ω que pertence àσ-álgebra M é dito mensurável. Finalmente, queremosassociar a cada evento do experimento uma probabili-dade de ocorrência. Uma medida de probabilidade é umafunção não negativa µ : M → [0, ∞) satisfazendo:

• µ(Ω) = 1;

• se (Ai)∞i=1 é uma família enumerável de elementos

disjuntos de M então

µ

∞

i=1

Ai

=

∞

∑i=1

µ(Ai) .

A tripla (Ω,M, µ) é um espaço de probabilidade e é ummodelo matemático conveniente para o tratamento demuitos sistemas probabilísticos. Voltando ao exemplodo arremesso de um dado, é fácil verificar que pondo

3

1 Algumas definições em probabilidade

1.1 Espaço de probabilidade

53Matemática Universitária nº46

Artigo

µ(i) = 1/6, para i = 1, . . . , 6, fica determinada umamedida de probabilidade.

1.2 Variavel Aleatoria

Num espaço de probabilidade (Ω,M, µ), conseguimosassociar uma probabilidade a qualquer evento mensu-rável. Em muitos modelos reais, gostaríamos de tra-tar matematicamente uma função real X : Ω → R,definida no espaço amostral Ω, isto é, calcular a pro-babilidade de eventos do tipo ω ∈ Ω : X(ω) < a,ω ∈ Ω : X(ω) ≤ a, ω ∈ Ω : X(ω) > a,ω ∈ Ω : X(ω) ≥ a etc. Para que isto seja possível,tais eventos precisam ser mensuráveis. Esta é a motiva-ção natural para o conceito de variável aleatória.

Seja um espaço de probabilidade (Ω,M, µ). Umafunção real X : Ω → R é dita mensurável se, para todoa ∈ R, ω ∈ Ω : X(ω) < a é mensurável, isto é, per-tence a M. Em probabilidade, uma função mensurávelé denominada uma variável aleatória. A partir desta defi-nição sobre R, estenderemos o conceito de variável ale-atória de forma natural sobre C, Cn e sobre o conjuntode matrizes complexas.

Uma função X : Ω → C é uma variável aleatória(complexa) se suas partes real e imaginária são variá-veis aleatórias (reais). Uma aplicação X : Ω → Cn

de coordenadas X = (X1, . . . , Xn) é um vetor aleató-rio se Xi : Ω → C é uma variável aleatória parai = 1, . . . , n. Denotaremos por Mm×n(C) o conjunto dasmatrizes m× n de números complexos. Uma aplicaçãoX = [Xij] : Ω → Mm×n(C) é uma matriz aleatória seXij : Ω → C é uma variável aleatória para i = 1, . . . , m ej = 1, . . . , n. Será conveniente identificar Cm com o con-junto Mm×1(C) das matrizes-coluna de números com-plexos. Neste sentido, vetores aleatórios são matrizes-coluna aleatórias.

1.3 Valor esperado

Dado um espaço de probabilidade (Ω,M, µ), um passocrucial num curso de teoria de medida é a definição daintegral

Xdµ de uma função mensurável X : Ω → R

relativa à medida µ. Não exporemos este procedimentoaqui por não ser necessário. O leitor interessado é con-vidado a consultar a ótima e rápida exposição em [3],

capítulos 4 e 5. Os fatos essenciais usados são as pro-priedades esperadas de linearidade e positividade daintegral:

•(aX + bY)dµ = a

Xdµ + b

Ydµ;

•

Xdµ ≥ 0 se X ≥ 0.

Quando|X|dµ < ∞, a função mensurável X : Ω → R

é dita integrável.Em probabilidade, o valor esperado E(X) da variável

aleatória X : Ω → R integrável é simplesmente a inte-gral de X: E(X) =

Xdµ.

Uma variável aleatória complexa X = Y + iZ : Ω →C é integrável se as partes real e imaginária são integrá-veis. O valor esperado de X é E(X) = E(Y) + iE(Z).Um vetor aleatório X : Ω → Cn é integrável se suascomponentes X1, . . . , Xn são variáveis aleatórias inte-gráveis, e o valor esperado de X é o vetor E(X) =(E(X1), . . . , E(Xn)). De forma análoga, dada uma ma-triz aleatória X = [Xij] : Ω → Mm×n(C) tal que suascomponentes Xij são variáveis aleatórias integráveis,define-se o valor esperado de X como a matriz com-plexa E(X) = [E(Xij)]. Daqui pra frente, sempre con-sideraremos variáveis aleatórias integráveis.

As propriedades de linearidade da integral se tra-duzem imediatamente para o operador valor esperadoE. Dadas as matrizes aleatórias X e Y, para quaisquerα, β ∈ C e quaisquer matrizes complexas A e B, temosque αX, βY, AX e XB definem matrizes aleatórias e tem-se E(αX + βY) = αE(X) + βE(Y), E(AX) = AE(X) eE(XB) = E(X)B, desde que as multiplicações de matri-zes façam sentido.

1.4 Media e covariancia

Como µ(Ω) = 1, para uma variável aleatória X tem-se E(X) =

Xdµ

µ(Ω) e, por isso, frequentemente dizemosque E(X) é a média da variável aleatória X. Analo-gamente diz-se que E(X) é a média do vetor aleató-rio X. Se X é uma variável aleatória real, o númeroVar(X) := E(X − E(X))2 é chamado variância de X.Para motivarmos uma interpretação melhor da variân-cia, informamos que f , g :=

f gdµ = E( f g) define

um produto interno no espaço vetorial das variáveisaleatórias reais f : Ω → R quadrado-integráveis, i.e.,tais que

| f |2dµ < ∞. Por sua vez, a norma induzida

4

1.3 Valor esperado

1.2 Variável Aleatória

1.4 Média e covariância

54 Matemática Universitária nº46

Artigo

por este produto interno é || f ||2 := f , f 12 = E(| f |2) 1

2 .Com esta noção de distância, vemos que a variânciaVar(X) = ||X − E(X)||22 mede o erro quadrático da va-riável aleatória em relação à sua média. Neste sentido,a variância dá uma interpretação de quanto a variávelaleatória se dispersa em relação à sua média. A gene-ralização para uma variável aleatória complexa é o nú-mero

Var(X) := E(X− E(X)) · (X− E(X))∗ ,

também chamado variância de X, em que denotamospor a∗ o complexo-conjugado do número a ∈ C. Taldefinição decorre do fato que, agora, o produto internoem questão é f , g :=

f g∗dµ = E( f g∗) no espaço

vetorial das variáveis aleatórias complexas f : Ω → C

quadrado-integráveis.O termo “valor esperado” ou “média” pode ser enga-

nador. Ele não significa, em geral, que a variável aleató-ria assume tipicamente este valor. Pode ser convenienteinterpretar a média de uma variável aleatória como,aproximadamente, a média de valores da variável alea-tória tomados sobre um número muito grande de ex-perimentos independentes. Aqui, o termo “indepen-dente” tem um significado preciso em probabilidade,mas que não abordaremos. O teorema que justifica talinterpretação da média de uma variável aleatória é cha-mado Lei dos Grandes Números, estudado em cursosde probabilidade. É importante citar que, em muitassituações em probabilidade e estatística, o mais impor-tante, ou pelo menos suficiente, é ter conhecimento damédia e da variância de uma variável aleatória.

A variância é uma medida de quanto a variável ale-atória se dispersa de sua própria média. Um conceitoque generaliza essa ideia para duas variáveis aleatóriasé a covariância. Dadas duas variáveis aleatórias X e Ydefinimos a covariância de X e Y como o número

Cov(X, Y) := E(X− E(X)) · (Y− E(Y))∗ .

A covariância entre duas variáveis aleatórias é uma me-dida de quanto as duas variáveis aleatórias variam con-juntamente, isto é, quanto a variação de uma em re-lação a sua média implica na variação da outra emrelação a sua média. Se Cov(X, Y) = 0 então dize-mos que X e Y são descorrelacionadas. Notemos queVar(X) = Cov(X, X).

Podemos agora generalizar tais conceitos para umvetor aleatório. A matriz complexa

Cov(X) := E(X− E(X)) · (X− E(X))∗

chama-se covariância do vetor aleatório X, em que X∗

é a transposta-conjugada de X. Futuramente, indicare-mos por AT a matriz transposta da matriz A. Notemosque se X = (X1, . . . , Xn) então Cov(X) = [Cov(Xi, Xj)],isto é, as entradas da matriz de covariância de um vetoraleatório são as covariâncias cruzadas de suas compo-nentes. Desta maneira, a matriz de covariância de umvetor aleatório carrega, de alguma forma, informaçãosobre como as componentes do vetor variam conjunta-mente.

Define-se a covariância-cruzada entre dois vetores ale-atórios X e Y como a matriz complexa

Cov(X, Y) := E(X− E(X)) · (Y− E(Y))∗ .

Novamente, se X = (X1, . . . , Xn) e Y = (Y1, . . . , Yn)então temos Cov(X, Y) = [Cov(Xi, Yj)]. Assim, a ma-triz covariância-cruzada Cov(X, Y) coleta as covariân-cias entre as componentes de X e Y. Se Cov(X, Y) = 0diz-se que X e Y são descorrelacionados, o que ocorrequando as componentes de X e de Y são descorrelaci-onadas.

Estamos prontos para a introdução de um poucode álgebra linear dentro desse contexto. Uma matrizcomplexa S ∈ Mm×m(C) é chamada de hermitiana seS = S∗. E é não negativa se, para qualquer vetorx ∈ Mm×1(C) vale x∗Sx ≥ 0. O lema seguinte usaum fato elementar de álgebra linear: se R é uma matrizcomplexa então RR∗ é matriz hermitiana não negativa(veja [5], pg 172).

Lema 1. A matriz de covariância S = E(X− E(X)) · (X−E(X))∗ de um vetor aleatório X é hermitiana não negativa.

Demonstração. Considere o vetor aleatórioX ∈ Mm×1(C). Tem-se S hermitiana, pois S∗ =[E(X − E(X)) · (X − E(X))∗]∗ = E[(X − E(X)) ·(X− E(X))∗]∗ = E(X− E(X))∗∗ · (X− E(X))∗ = S.Suponha E(X) = 0, para facilitar. Dado um vetor x ∈Mm×1(C), tem-se x∗Sx = x∗E(XX∗)x = E(x∗XX∗x).Tem-se Y = x∗XX∗x ≥ 0, pois XX∗ é hermitianonão negativo. Assim, Y é uma variável aleatória não

5

55Matemática Universitária nº46

Artigo

negativa, logo E(Y) =

Ydµ ≥ 0, já que a integralé um operador não negativo. Segue que x∗Sx ≥ 0.Portanto, S é hermitiana não negativa. O caso E(X) = 0se resolve da mesma maneira para o vetor aleatórioY = X− E(X).

No início do artigo, apresentamos o problema das se-noides com a promessa de formalizá-lo. Agora que te-mos os ingredientes necessários da teoria de probabili-dade, cumpriremos a promessa.

Temos um espaço de probabilidade com espaçoamostral Ω. No sistema de m equações

Xk =p

∑j=1

γjeiθjk + Wk, k = 1, . . . , m , (1)

X1, . . . , Xm, W1, . . . , Wm e γ1, . . . , γp são variáveis ale-atórias complexas definidas em Ω e θ1, . . . , θp são nú-meros reais. Nosso objetivo é determinar os númerosθ1, . . . , θp usando apenas as propriedades estatísticasdas variáveis aleatórias X1, . . . , Xm medidas no recep-tor, isto é, suas médias, variâncias e covariâncias.

O sistema acima pode ser reescrito numa equaçãomatricial mais compacta. Em primeiro lugar, definimosos vetores

X = [X1 · · ·Xm]T e W = [W1 · · ·Wm]T .

Defina a aplicação a : R → Mm×1(C) de componen-tes a = [a1 · · · am]T tal que ak(θ) = ei(k−1)θ , ∀θ ∈ R ek = 1, . . . , m. Em outras palavras

a(θ) =1 eiθ ei2θ · · · ei(m−1)θ

T, ∀θ ∈ R .

Para facilitar, definiremos o vetor

θ := (θ1, . . . , θp)

cujas componentes são as frequências das senoides a se-rem determinadas. A partir da aplicação a definimos amatriz

Aθ

:= [ak(θj)] ∈ Mm×p(C) ,

ou seja, Aθ

é a matriz cujas colunas são os vetoresa(θ1), . . . , a(θp). Note que, como explicitado, essa ma-triz é função do vetor θ a determinar.

Finalmente, defina também o vetor F = [F1 · · · Fp]T

com Fj = γjeiθj , para j = 1, . . . , p, isto é:

F =γ1eiθ1 · · · γpeiθp

T.

Em resumo, com esta notação ficam definidos os ve-tores aleatórios

X : Ω → Mm×1(C),W : Ω → Mm×1(C),F : Ω → Mp×1(C),

ficando a cargo do leitor verificar que o sistema deequações (1) é equivalente à equação matricial

X = Aθ

F + W . (2)

A seguir, listaremos as hipóteses a serem satisfeitaspara que possamos solucionar o problema das senoi-des.

Em primeiro lugar, quanto às características estatísti-cas do ruído, supomos que

E1: W tem média zero e covariância Cov(W) = σ2Ipara algum σ > 0.

Acima, I é a matriz identidade de Mm×m(C). Neste ar-tigo, chamaremos de ruído branco gaussiano o vetor ale-atório satisfazendo a hipótese acima. Desta hipótese re-sulta que todas as componentes de W têm média zeroe que componentes diferentes de W são descorrelaci-onadas, o que nos informa que elas variam de forma“independente”. Além disso, todas as componentes deW têm variância σ2, isto é, as componentes de W estãoigualmente dispersas de sua média (zero). Surpreen-dentemente, vários contextos reais podem ser modela-dos por um ruído gaussiano de forma bem sucedida.Isto advém de um resultado profundo em probabili-dade chamado Teorema do Limite Central.

Também supomos que

E2: W e F são descorrelacionados.

Ou seja, Cov(W, F) = 0. A interpretação desta hipó-tese é que as variações aleatórias do ruído W e do sinaltransmitido F não têm correlação, o que não é uma hi-pótese muito forte em vários contextos.

A seguir, duas condições técnicas serão necessárias.A primeira é

6

2 O problema

56 Matemática Universitária nº46

Artigo

H1: m > p.

Apesar desta hipótese ser tecnicamente essencial no de-correr da solução, ela possui uma interpretação: o nú-mero m de unidades de tempo do sinal medido pelo re-ceptor é suficientemente grande em relação ao númerop de parâmetros a serem determinados. Do ponto devista prático, isto é, de certa forma, intuitivo. A se-gunda hipótese técnica é sobre a característica estatís-tica do sinal transmitido F:

H2: Cov(F) é inversível.

Um cálculo simples mostra que a matriz Cov(F), a co-variância de F, é a matriz cuja entrada de linha k e co-luna l é Cov(γk, γl)ei(θk−θl). Desta forma, se definimosγ :=

γ1, . . . , γp

T , então Cov(F) depende da matrizCov(γ) e podemos pensar na hipótese H2 como umacondição sobre γ, isto é, sobre as amplitudades das se-noides. A hipótese H2 é forte demais em alguns con-textos. Neste caso, o algoritmo MUSIC não pode solu-cionar o problema das senoides e outros métodos sãonecessários.

Passamos agora às condições sobre a aplicação a. No-temos que a : R → Mm×1(C) é uma aplicação contínua.Além disso, a é periódica de período 2π. Para determi-narmos as frequências θ1 . . . , θp, usando a aplicação a,necessitaremos de sua injetividade, já que a estratégiaserá explorar os vetores a(θ1), . . . , a(θp) colunas da ma-triz A

θ. Além disso, não é esperado que consigamos

determinar duas frequências iguais. Por esta razão, de-notando por U o intervalo aberto (0, 2π), supomos que

A1: θ1, . . . , θp são dois a dois distintos e estão em U.

Convém notar que a hipótese A1 não é uma restrição,pois, módulo 2π, sempre podemos reescrever as senoi-des de modo a satisfazer tal condição. Notemos que

A2: a restrição a|U é injetiva.

Para facilitar, de agora em diante denotaremos a matrizA

θ

simplesmente por A. Pediremos que

A3: posto(A) = p.

Note que esta condição é automaticamente satisfeita segarantirmos as hipóteses H1 e A1. De fato, se θ1, . . . , θp

são distintos e contidos no intervalo (0, 2π) temos que

A =

1 1 · · · 1eiθ1 eiθ2 · · · eiθp

ei2θ1 ei2θ2 · · · ei2θp

......

......

ei(m−1)θ1 ei(m−1)θ2 · · · ei(m−1)θp

é uma matriz de Vandermond. Logo as colunas de Asão linearmente independentes. Como m > p, isto é,o número de linhas é maior que o número de colunas,segue que posto(A) = p. Finalmente, notemos que asatisfaz a seguinte condição geométrica:

A4: a(U) ∩ [a(θ1), . . . , a(θp)] = a(θ1), . . . , a(θp),

em que [a(θ1), . . . , a(θp)] denota o subespaço geradopelos vetores a(θ1), . . . , a(θp). A utilidade de A4 re-side no fato de que, dado θ ∈ U, basta que a(θ) es-teja em [a(θ1), . . . , a(θp)] para garantir que θ é umadas frequências θ1, . . . , θp procuradas. Para mostrar A4,basta proceder de forma análoga à prova de A3. Defato, se θ ∈ U − θ1, . . . , θp então a matriz

1 1 1 · · · 1eiθ eiθ1 eiθ2 · · · eiθp

ei2θ ei2θ1 ei2θ2 · · · ei2θp

......

......

...ei(m−1)θ ei(m−1)θ1 ei(m−1)θ2 · · · ei(m−1)θp

também é matriz de Vandermond. Logo suas colunasa(θ), a(θ1), . . . , a(θp) são linearmente independentes esegue que a(θ) /∈ [a(θ1), . . . , a(θp)].

Sob estas condições, o problema das senoides con-siste em determinar os p números θ1, . . . , θp tendo co-nhecimento apenas da matriz de covariância de X e daaplicação a|U, tal que os vetores aleatórios em questãosatisfaçam a equação linear (2).

Por completude, afirmamos que o algoritmo MUSICpode ser utilizado num contexto um pouco mais geraldo que o problema das senoides. De fato, basta que omodelo seja a equação (2) e que as hipóteses listadassejam satisfeitas. Entretanto, o exemplo das senoides ébastante ilustrativo e já possui vasta aplicação.

7

57Matemática Universitária nº46

Artigo

Iremos agora solucionar o problema proposto. Isto seráfeito construindo-se uma função f : U → [0, ∞) realnão negativa no aberto U, cujos p zeros são exatamenteθ1, . . . , θp ∈ R. Além disso, f será construída usando-sesomente a matriz de covariância S := E(X − E(X)) ·(X − E(X))∗ do vetor aleatório X e a aplicação a res-trita a U.

A construção da função f será geométrica. Mais pre-cisamente, iremos construir a partir de S um subes-paço vetorial Vσ2 de Cm que é complemento ortogo-nal do subespaço [a(θ1), . . . , a(θp)] gerado pelos vetoresa(θ1), . . . , a(θp). Este é o conteúdo da Proposição 4. Aconstrução de f segue dessa proposição de forma ime-diata projetando-se, para cada θ ∈ U, o valor a(θ) sobreVσ2 e observando-se que θ ∈ U é um dos parâmetrosθ1, . . . , θp se e somente se a(θ) está em [a(θ1), . . . , a(θp)].Este é o conteúdo do teorema final.

Resta a pergunta: como determinamos o espaço ve-torial Vσ2? O que faremos será mostrar que o parâme-tro σ2 associado ao ruído é o menor autovalor da matrizS e Vσ2 será o autoespaço de S associado ao autovalorσ2 (lembre-se que sendo S hermitiana, seus autovalo-res são reais e portanto podem ser ordenados). Este é oconteúdo das Proposições 1, 2, e 3.

Além da matriz de covariância S, definimos também

P := E(F− E(F)) · (F− E(F))∗ .

A seguir, dois resultados auxiliares importantes são ob-tidos.

Lema 2. posto(AP) = posto(APA∗) = p.

Demonstração. Para tanto precisaremos de um resul-tado básico de álgebra linear (veja [5], exercício 8.41, p.104): dadas duas matrizes R e T, tem-se posto(RT) ≤posto(R) e posto(RT) ≤ posto(T), desde que a multi-plicação faça sentido. Temos, portanto, que posto(A) ≥posto(AP). Por hipótese, P é inversível, de modo que[AP]P−1 = A e, portanto, posto(A) ≤ posto(AP). Con-cluímos que posto(AP) = posto(A) = p. De formaanáloga, tem-se a desigualdade posto(A[PA∗]) ≤posto(A) = p. AP é matriz m× p com posto(AP) = p,portanto AP tem inversa à esquerda, que denotaremos

por B. Assim, B[APA∗] = A∗ e, portanto,

posto(APA∗) ≥ posto(A∗) = posto(A) = p .

Concluímos que posto(APA∗) = p.

A equação (2) que rege o problema é escrita direta-mente em termos dos vetores aleatórios X, F e W. Umaequação análoga mais útil em termos das covariânciasde tais vetores aleatórios pode ser obtida. Note que acondição m > p implica que a matriz APA∗ é singu-lar, pois APA∗ é matriz quadrada m × m de posto p.Essa condição permitirá relacionar os autovalores dasmatrizes S e APA∗. Este é o conteúdo da proposiçãoseguinte.

Proposição 1. Tem-se

S = APA∗ + σ2I . (3)

O número σ2 é autovalor de S e as matrizes S e APA∗ pos-suem os mesmos autoespaços. Finalmente, se λ é autovalorde S e µ é autovalor de APA∗ cujos autoespaços são iguaisentão

λ = µ + σ2 . (4)

Demonstração. Para facilitar chame X := X − E(X) eF := F− E(F). Da equação (2) tem-se

X = AF + W−AE(F) = AF + W ,

já que W tem média zero. Usando a linearidade de E,tem-se

S = EX · X∗

= E(AF + W)(AF + W)∗

= AE(F · F∗)A∗ + AE(FW∗) +

E(WF∗)A∗ + E(WW∗).

Usando-se que E(FW∗) = E(WF∗) = 0, P = E(F ·F∗) e E(WW∗) = σ2I, segue que S = APA∗ + σ2I.Sendo APA∗ singular, det(APA∗) = 0 e, usando-se (3),conclui-se que det(S− σ2I) = 0. Logo, σ2 é autovalorde S. A relação (4) é obtida sem maiores dificuldadesusando-se a equação (3). De fato, denote por Vλ o auto-espaço de S associado a λ e Uµ o autoespaço de APA∗

associado a µ. Então

e ∈ Vλ ⇔ (APA∗ + σ2I)e = λe ⇔

APA∗e = (λ− σ2)e ⇔ e ∈ Uλ−σ2 .

8

3 A solução

58 Matemática Universitária nº46

Artigo

As equações (3) e (4) são equações análogas à equa-ção (2) em outros termos. A primeira é dada em termosdas matrizes de covariância e a segunda em termos dosautovalores das matrizes de covariância e da matriz au-xiliar APA∗.

Note que mesmo que saibamos a partir da Proposição1 que o parâmetro σ2 associado ao ruído é um autova-lor de S, não conseguimos determiná-lo (S pode ter mautovalores distintos). Isto é o que faremos agora.

Proposição 2 (Estimativa do ruído). O número σ2 é o me-nor autovalor da matriz de covariância S.

Demonstração. Usaremos basicamente o fato de que Se APA∗ são matrizes hermitianas não negativas e quepossuem seus autovalores relacionados por (4), con-forme a Proposição 1. S é hermitiana não negativa poisé matriz de covariância (Lema 1). Denotando nova-mente F := F− E(F), temos

APA∗ = AE(F · F∗)A∗ = E(AF · F∗A∗)

= E[AF(AF)∗] = E(RR∗),

em que R := AF ∈ Mm×1(C). Como E(R) = AE(F) =0, pois E(F) = 0, segue que APA∗ é a matriz de co-variância do vetor aleatório R. Sendo assim, do Lema 1segue que APA∗ é hermitiana não negativa. Portanto, Se APA∗ possuem todos autovalores reais não negativosrelacionados por (4). Denote os autovalores de S porλi(S), para i = 1, . . . , k, e denote os respectivos autova-lores de APA∗, segundo a equação (4), por λi(APA∗),para i = 1, . . . , k. Isto significa que λi(S) = λi(APA∗)+σ2, para i = 1, . . . , k. Como APA∗ é singular, possui au-tovalor nulo; segue que mini=1,...,k λi(APA∗) = 0. Se-gue da equação (4) que

mini=1,...,k

λi(S) = mini=1,...,k

λi(APA∗) + σ2 = σ2 .

Observe que também não é necessário fornecer o nú-mero de parâmetros p de interesse. De fato, este númeroé obtido pela proposição seguinte.

Proposição 3 (Estimativa do número de parâmetros).Se n é a multiplicidade algébrica (igual à multiplicidade geo-métrica) do menor autovalor σ2 da matriz de covariância S,então

p = m− n .

Demonstração. Dada uma matriz R ∈ Ma×b(C), a se-guir denotaremos por N(R) = x ∈ Mb×1(C) : Rx =0 o núcleo da matriz R. A mesma notação é usadapara transformações lineares. Usaremos que (veja [5],pg 287), para uma matriz hermitiana R, a multiplici-dade algébrica de todo autovalor λ é igual à sua mul-tiplicidade geométrica (dimensão do autoespaço Vλ =N(R− λI) formado pelos autovetores associados a λ).

Pela Proposição 1, S e APA∗ possuem os mesmos au-toespaços, de modo que, se λ(S) é autovalor de S eλ(APA∗) é autovalor de APA∗, os quais se correspon-dem pela equação (4), então λ(S) e λ(APA∗) possuem amesma multiplicidade geométrica. APA∗ possui o au-tovalor nulo λ0 = 0 (por ser singular) que correspondeao menor autovalor σ2 de S (Proposição 2), segundo aequação (4). Sendo S hermitiana, a multiplicidade geo-métrica do autovalor σ2 é n. Logo a multiplicidade geo-métrica de λ0 é n. Mas a multiplicidade geométricade λ0 = 0 é dim N(APA∗ − 0I) = dim N(APA∗) =dim Mm(C)− posto(APA∗) = m− p, pelo Teorema doNúcleo e Imagem e o Lema 2. Segue, portanto, quen = m− p.

A seguir, por uma questão de gosto, trabalhare-mos com transformações lineares, ao invés de matri-zes. Neste sentido, as matrizes A ∈ Mm×p(C), P ∈Mp×p(C) e S ∈ Mm×m(C) são, respectivamente, matri-zes das transformações lineares

A : Mp×1(C)→ Mm×1(C),P : Mp×1(C)→ Mp×1(C),S : Mm×1(C)→ Mm×1(C),

relativamente às respectivas bases canônicas. Comousamos as bases canônicas, os autovalores e autoveto-res da matriz S são os mesmos da transformação linearS e a transformação linear S é hermitiana não negativa,pois a matriz S o é (a base canônica é ortonormal). Omesmo pode ser dito sobre a matriz APA∗ e a transfor-mação APA∗.

O conjunto a(θ1), . . . , a(θp) ⊂ Mm×1(C) das colu-nas de A está contido na imagem Im(A) da transfor-mação A. Mas o conjunto a(θ1), . . . , a(θp) é linear-mente independente, pois é formado pelas colunas deA e posto(A) = p. Portanto, a(θ1), . . . , a(θp) é, defato, uma base de Im(A).

9

59Matemática Universitária nº46

Artigo

Pelo Teorema Espectral, existe uma base ortonormalde Mm×1(C) formada por autovetores de S. Mais pre-cisamente, Mm×1(C) pode ser decomposto na soma di-reta dos autoespaços de S, isto é,

Mm×1(C) = Vλ1 ⊕ · · · ⊕Vλk ,

em que λ1, . . . , λk são os autovalores de S e, para cadai = 1, . . . , k, Vλi é o autoespaço de S associado ao au-tovalor λi. O menor autovalor σ2 de S (Proposição 2)tem multiplicidade algébrica n, que é igual à sua mul-tiplicidade geométrica, pois S é hermitiana. Portanto,dim Vσ2 = n. Temos que Im(A) e Vσ2 são subespa-ços de Mm×1(C). Mostraremos na Proposição 4 que,na verdade, esses dois subespaços são complementosortogonais um do outro, isto é, Im(A)⊥ = Vσ2 . Essapropriedade geométrica será fundamental na determi-nação dos parâmetros de interesse θ1, . . . , θp.

Proposição 4. Im(A)⊥ = Vσ2 .

Demonstração. Usaremos um fato elementar de álgebralinear: se R : E → F é transformação linear entredois espaços vetoriais E e F com produto interno entãoN(R∗) = Im(R)⊥ (veja [5], pg 144).

Em primeiro lugar, a transformação linear AP :Mp×1(C) → Mm×1(C) tem posto(AP) = p, pelo Lema2, logo AP é injetiva e portanto tem inversa à esquerdaB : Mm×1(C) → Mp×1(C). Seja u1, . . . , un uma baseortonormal de Vσ2 . Os autovetores u1, . . . , un de S as-sociados ao autovalor σ2 são os autovetores de APA∗

associados ao autovalor nulo, segundo a Proposição 1.Segue que, para i = 1, . . . , n, tem-se APA∗ · ui = 0 ⇒BAPA∗ · ui = 0 ⇒ A∗ · ui = 0. Essa afirmação significaque u1, . . . , un ⊂ N(A∗) = Im(A)⊥. Concluímos queVσ2 ⊂ Im(A)⊥, logo Vσ2 é subespaço de Im(A)⊥. Mas

dim Im(A)⊥ = dim Mm×1(C)− posto(A) =

= m− p = n = dim Vσ2 ,

pela Proposição 3. Portanto, Vσ2 = Im(A)⊥.

Estamos preparados para determinar os parâmetrosde interesse através da Proposição 4. Considere no-vamente uma base ortonormal u1, . . . , un de Vσ2 . See1, . . . , en é a base canônica de Mn×1(C), existe trans-formação linear T : Mn×1(C) → Mm×1(C) tal que

Te1 = u1, . . . , Ten = un. Relativamente às bases canô-nicas de Mn×1(C) e Mm×1(C), a matriz T ∈ Mm×n(C)dessa transformação linear é a matriz cujas colunas sãoos vetores u1, . . . , un ∈ Mm×1(C). Note que T é inje-tiva e Im(T) = Vσ2 .

Teorema. Os parâmetros θ1, . . . , θp são os p zeros da fun-ção f : U → [0, ∞) dada por

f (θ) = T∗ · a(θ)2 .

Em outros termos, seja g = 1/ f : U − f−1(0) → (0, ∞)dada por

g(θ) =1

T∗ · a(θ)2 . (5)

Para todo θ ∈ U, θ ∈ θ1, . . . , θp se e somente se

limx→θ

g(x) = ∞ .

Demonstração. Note que, dado θ ∈ U, f (θ) = 0 see somente se a(θ) ∈ Im(A). De fato, f (θ) = 0 ⇔T∗ · a(θ) = 0 ⇔ T∗ · a(θ) = 0 ⇔ a(θ) ∈ N(T∗) =Im(T)⊥ = V⊥

σ2 ⇔ a(θ) ∈ Im(A), pela Proposição 4.

Tem-se trivialmente que θ1, . . . , θp ⊂ f−1(0), poisa(θi) ∈ Im(A). Reciprocamente, seja θ ∈ U comf (θ) = 0, então a(θ) ∈ Im(A). Por hipótese temosθ1, . . . , θp ⊂ U e a(U) ∩ Im(A) = a(θ1), . . . , a(θp),já que o espaço gerado pelos vetores a(θ1), . . . , a(θp)é Im(A). Portanto, a(θ) ∈ a(θ1), . . . , a(θp), ou seja,a(θ) = a(θi) para algum i ∈ 1, . . . , p. Sendo a|U inje-tiva, segue que θ = θi, isto é, θ ∈ θ1, . . . , θp. Concluí-mos que f−1(0) = θ1, . . . , θp.

Analisemos a função g. Se a aplicação a é contínuaem U então a função f é contínua por ser compostade funções contínuas e, também, g é contínua. Se i ∈1, . . . , p, tem-se f (θi) = 0. Sabe-se que f−1(0) é dis-creto, por ser um subconjunto finito de R. Logo, existeuma bola aberta B(θi, ri) ⊂ U tal que f é estritamentepositiva em B(θi, ri) − θi. Usando isto e a continui-dade de f em θi, segue que limx→θi f (x) = f (θi) = 0⇒ limx→θi g(x) = limx→θi

1f (x) = ∞. Reciprocamente,

seja θ ∈ U − θ1, . . . , θp = U − f−1(0), então g é con-tínua em θ, logo limx→θ g(x) = g(θ) < ∞. Concluímosque para todo θ ∈ U, θ ∈ θ1, . . . , θp se e somentelimx→θ g(x) = ∞.

10

60 Matemática Universitária nº46

Artigo

Analisando o teorema anterior, o leitor pode constatarque temos em mãos uma função ( f ou g) que permitedeterminar as frequências procuradas θ1, . . . , θp das se-noides. Repare que mantivemos a promessa de usarapenas a matriz de covariância S do vetor aleatório Xe a aplicação a|U: f e g são construídas a partir da apli-cação a|U e da matriz T (recorde a contrução de T: ascolunas de T são uma base ortonormal de Vσ2 , o autoes-paço da aplicação S associado ao seu menor autovalorσ2). A função g é mais utilizada em aplicações com-putacionais pois “picos” são mais robustos a erros decálculo do que zeros.

Este artigo descreveu um problema de probabilidadecuja solução foi construída usando-se essencialmenteresultados de álgebra linear. Quanto ao âmbito de apli-cações, um detalhe é crucial: raramente se conheceo valor exato da matriz de covariância S = E(X −E(X))(X− E(X))∗ do vetor aleatório X. Isto se deve aofato de que a medição experimental permite determinaruma coleção finita X(ωi) : ωi ∈ Ω, i = 1, . . . , K devalores do vetor aleatório X, em que cada valor X(ωi)(chamado de realização do vetor aleatório X) é obtidono i-ésimo experimento. Neste caso, conforme men-cionado brevemente na introdução do artigo, define-se a matriz covariância estimada S := 1

K ∑Ki=1(X(ωi) −

E(X))(X(ωi) − E(X))∗, em que E(X) := 1K ∑K

i=1 X(ωi)é a média estimada, e pergunta-se em que condições esta(ou outras estimativas) pode substituir a “verdadeira”matriz de covariância S. Essa é a motivação para a teo-ria de estimação (veja por exemplo [4], capítulo 2) quepertence à estatística matemática. Apesar de não mos-trarmos aqui, uma consequência importante de usar-mos a matriz S ao invés de S é que o número p de pa-râmetros de interesse precisa também ser estimado (oque era desnecessário com o conhecimento de S, videProposição 3).

Esperamos que o estudante de matemática possa terapreciado com este exemplo uma aplicação elegante daálgebra linear em outros contextos da matemática e damatemática aplicada.

Agradecimentos. Agradeço ao parecerista e aos edi-tores pelas inúmeras sugestões para a melhoria deste

texto, e também ao CNPq, que me financiou durante omestrado.

[1] STELZER, A.; STELZER, M. Resolution enhance-ment with model-based frequency estimation al-gorithms in radar signal processing. Surface Sen-sing Technologies and Applications, v. 4, n. 3, p. 241–261, 2003.

[2] THOMPSON, P.; NANNINI, M.; SCHEIBER, R. Tar-get separation in SAR image with the MUSIC algo-rithm. In: IEEE International Symposium in Geos-cience and Remote Sensing, Barcelona, 2007. Pro-ceedings. Barcelona, IEEE, 2007. p. 468–471.

[3] BARTLE, R. G. Elements of integration. New York:John Wiley, 1966.

[4] BOROVKOV, A. A. Mathematical statistics. Amster-dam: Gordon and Breach, 1998.

[5] LIMA, E. L.. Álgebra linear. 7 ed. Rio de Janeiro:IMPA, 2006. (Coleção Matemática Universitária)

[6] PEEBLES, P. Z. Probability, random variables, and ran-dom signal principles. 2 ed. New York: McGraw-Hill, 1987.

[7] SCHMIDT, R. O. Multiple emitter location and sig-nal parameter estimation. IEEE Transactions on An-tennas and Propagation, AP-34, n. 10, p. 276–280,1986.

Philip Thompsonphilipthomp@gmail.com

11

4 conclusão

Referências