Analyse avec des données complexes, quelque chose de différent?

31

Disons par exemple que vous faites un modèle linéaire, mais les données sont complexes. $y$

$y = x \beta + \epsilon$

Mon ensemble de données est complexe, comme dans tous les nombres en sont de la forme . Y a-t-il quelque chose de différent sur le plan de la procédure lorsque vous travaillez avec de telles données? $y$ $(a + bi)$

Je demande parce que vous finirez par obtenir des matrices de covariance complexes et des statistiques de test qui sont complexes.

Avez-vous besoin d'utiliser un transposé conjugué au lieu de transposer lorsque vous faites des moindres carrés? une covariance de valeur complexe est-elle significative?

regression anova data-transformation complex-numbers bill_e
la source

3

Considérez le nombre complexe comme deux variables distinctes, et supprimez ainsi i de toutes vos équations. Sinon, ce sera un cauchemar ...

sashkello

Des informations sur

ou

?

x

$x$

β

$\beta$

Stijn

3

@Sashkello Quel "cauchemar"? Les dimensions sont divisées par deux lorsque vous utilisez des nombres complexes, donc c'est sans doute une simplification. De plus, vous avez transformé un DV bivarié en DV univarié , ce qui est un énorme avantage. Peter Rabbit: oui, des transpositions conjuguées sont nécessaires. La matrice de covariance complexe est hermitienne positive-définie. Comme son homologue réel, il a toujours des valeurs propres réelles positives, ce qui répond à la question du sens.

whuber

2

@whuber Cela n'a aucun sens pour moi d'entrer dans des nombres complexes si le problème est comme indiqué. Il n'est pas plus simple de traiter des nombres complexes - sinon il n'y aurait pas du tout de question ici. Tout ne fonctionnera pas bien avec des nombres complexes et ce n'est pas un changement simple si vous ne savez pas ce que vous faites. Transformer ce problème dans l'espace réel est équivalent , et vous pouvez alors appliquer toute la variété des techniques statistiques sans vous inquiéter si cela fonctionne ou non dans un espace complexe.

sashkello le

1

@whuber Bonne réponse et belle explication. Je dirais que dès que vous passez la transformation de l'un à l'autre, ce n'est vraiment pas difficile ...

sashkello

40

Sommaire

La généralisation de la régression des moindres carrés aux variables à valeurs complexes est simple, consistant principalement à remplacer les transpositions matricielles par des transpositions conjuguées dans les formules matricielles habituelles. Une régression à valeurs complexes, cependant, correspond à une régression multiple multivariée compliquée dont la solution serait beaucoup plus difficile à obtenir en utilisant des méthodes standard (variables réelles). Ainsi, lorsque le modèle à valeurs complexes est significatif, l'utilisation d'une arithmétique complexe pour obtenir une solution est fortement recommandée. Cette réponse comprend également des suggestions de façons d'afficher les données et de présenter des tracés de diagnostic de l'ajustement.

Par souci de simplicité, discutons le cas de la régression ordinaire (univariée), qui peut être écrite

z_{j} = β_{0} + β_{1} w_{j} + ε_{j} .

$z_j = \beta_0 + \beta_1 w_j + \varepsilon_j.$

J'ai pris la liberté de nommer la variable indépendante et la variable dépendante , qui est conventionnelle (voir, par exemple, Lars Ahlfors, Complex Analysis ). Tout ce qui suit est simple à étendre au paramètre de régression multiple. $W$ $Z$

Interprétation

Ce modèle a une interprétation géométrique facilement visualisable: la multiplication par redimensionnera par le module de et le fera tourner autour de l'origine par l'argument de . Par la suite, l'ajout de traduit le résultat par ce montant. L'effet de est de "gigue" cette traduction un peu. Ainsi, régresser le sur le de cette manière est un effort pour comprendre la collection de points 2D $\beta_1$ $w_j$ $\beta_1$ $\beta_1$ $\beta_0$ $\varepsilon_j$ $z_j$ $w_j$ $(z_j)$ comme résultant d'une constellation de points 2D via une telle transformation, permettant une erreur dans le processus. Ceci est illustré ci-dessous avec la figure intitulée "Fit as a Transformation". $(w_j)$

Notez que le redimensionnement et la rotation ne sont pas n'importe quelle transformation linéaire du plan: ils excluent par exemple les transformations de biais. Ainsi , ce modèle n'est pas la même que celle d' une régression multiple avec deux variables quatre paramètres.

Les moindres carrés ordinaires

Pour connecter le cas complexe avec le cas réel, écrivons

pour les valeurs de la variable dépendante et $z_j = x_j + i y_j$

pour les valeurs de la variable indépendante. $w_j = u_j + i v_j$

De plus, pour les paramètres, écrivez

et . $\beta_0 = \gamma_0 + i \delta_0$ $\beta_1 = \gamma_1 +i \delta_1$

Chacun des nouveaux termes introduits est, bien sûr, réel, et est imaginaire tandis que indexe les données. $i^2 = -1$ $j=1, 2, \ldots, n$

OLS découvertes et qui minimisent la somme des carrés des écarts, $\hat\beta_0$ $\hat\beta_1$

\sum_{j = 1}^{n} | | z_{j} - ({\hat{β}}_{0} + {\hat{β}}_{1} w_{j}) | |^{2} = \sum_{j = 1}^{n} ({\bar{z}}_{j} - (\bar{{\hat{β}}_{0}} + \bar{{\hat{β}}_{1}} {\bar{w}}_{j})) (z_{j} - ({\hat{β}}_{0} + {\hat{β}}_{1} w_{j})) .

$\sum_{j=1}^n ||z_j - \left(\hat\beta_0 + \hat\beta_1 w_j\right)||^2 = \sum_{j=1}^n \left(\bar z_j - \left(\bar{\hat\beta_0} + \bar{\hat\beta_1} \bar w_j\right)\right) \left(z_j - \left(\hat\beta_0 + \hat\beta_1 w_j\right)\right).$

Formellement, cela est identique à la formulation matricielle habituelle: comparez-la à La seule différence que nous trouvons est que la transposée de la matrice de conception est remplacée par la transposée conjuguée . Par conséquent, la solution matricielle formelle est $\left(z - X\beta\right)'\left(z - X\beta\right).$ $X'$ $X^* = \bar X '$

\hat{β} = {(X^{*} X)}^{- 1} X^{*} z .

$\hat\beta = \left(X^*X\right)^{-1}X^* z.$

Dans le même temps, pour voir ce qui pourrait être accompli en transposant cela en un problème purement variable, nous pouvons écrire l'objectif OLS en termes de composants réels:

\sum_{j = 1}^{n} {(x_{j} - γ_{0} - γ_{1} u_{j} + δ_{1} v_{j})}^{2} + \sum_{j = 1}^{n} {(y_{j} - δ_{0} - δ_{1} u_{j} - γ_{1} v_{j})}^{2} .

$\sum_{j=1}^n \left(x_j-\gamma_0-\gamma_1u_j+\delta_1v_j\right)^2 + \sum_{j=1}^n\left(y_j-\delta_0-\delta_1u_j-\gamma_1v_j\right)^2.$

Cela représente évidemment deux régressions réelles liées : l'une régresse sur et , l'autre régresse sur et ; et nous exigeons que le coefficient pour soit le négatif du coefficient pour et que le coefficient pour égal au coefficient pour . De plus, parce que le total $x$ $u$ $v$ $y$ $u$ $v$ $v$ $x$ $u$ $y$ $u$ $x$ $v$ $y$ les carrés des résidus des deux régressions doivent être minimisés, il ne sera généralement pas le cas que l'un ou l'autre ensemble de coefficients donne la meilleure estimation pour ou seul. Ceci est confirmé dans l'exemple ci-dessous, qui réalise séparément les deux régressions réelles et compare leurs solutions à la régression complexe. $x$ $y$

Cette analyse montre que la réécriture de la régression complexe en termes de parties réelles (1) complique les formules, (2) obscurcit l'interprétation géométrique simple et (3) nécessiterait une régression multiple multivariée généralisée (avec des corrélations non triviales entre les variables ) résoudre. On peut faire mieux.

Exemple

À titre d'exemple, je prends une grille de valeurs aux points intégraux près de l'origine dans le plan complexe. Aux valeurs transformées s'ajoutent des erreurs ayant une distribution gaussienne bivariée: en particulier, les parties réelles et imaginaires des erreurs ne sont pas indépendantes. $w$ $w\beta$

Il est difficile de tracer le diagramme de dispersion habituel de pour les variables complexes, car il serait composé de points en quatre dimensions. Au lieu de cela, nous pouvons voir la matrice de nuage de points de leurs parties réelles et imaginaires. $(w_j, z_j)$

Matrice de nuage de points

Ignorez l'ajustement pour l'instant et regardez les quatre premières lignes et les quatre colonnes de gauche: elles affichent les données. La grille circulaire de est évidente dans le coin supérieur gauche; il compte points. Les diagrammes de dispersion des composantes de rapport aux composantes de montrent des corrélations claires. Trois d'entre eux ont des corrélations négatives; seuls (la partie imaginaire de ) et (la partie réelle de ) sont corrélés positivement. $w$ $81$ $w$ $z$ $y$ $z$ $u$ $w$

Pour ces données, la vraie valeur de est $\beta$ . Elle représente une extension deet une rotationsens antihoraire de 120 degréssuivi partraduction deunités vers la gauche etunités vershaut. Je calcule trois ajustements: la solution des moindres carrés complexes et deux solutions OLS pouretséparément, pour comparaison. $(-20 + 5i, -3/4 + 3/4\sqrt{3}i)$ $3/2$ $20$ $5$ $(x_j)$ $(y_j)$

Fit            Intercept          Slope(s)
True           -20    + 5 i       -0.75 + 1.30 i
Complex        -20.02 + 5.01 i    -0.83 + 1.38 i
Real only      -20.02             -0.75, -1.46
Imaginary only          5.01       1.30, -0.92

Il sera toujours le cas que l'ordonnée à l'origine uniquement en accord avec la partie réelle de l'ordonnée à l'origine complexe et l'ordonnée à l'origine uniquement en accord avec la partie imaginaire de l'ordonnée à l'origine complexe. Il est évident, cependant, que les pentes réelles et imaginaires uniquement ne correspondent ni aux coefficients de pente complexes ni entre elles, exactement comme prévu.

Examinons de plus près les résultats de l'ajustement complexe. Tout d'abord, un tracé des résidus nous donne une indication de leur distribution gaussienne bivariée. (La distribution sous-jacente a des écarts-types marginaux de et une corrélation de .) Ensuite, nous pouvons tracer les magnitudes des résidus (représentés par les tailles des symboles circulaires) et leurs arguments (représentés par des couleurs exactement comme dans le premier tracé) par rapport aux valeurs ajustées: ce tracé devrait ressembler à une distribution aléatoire de tailles et de couleurs, ce qu'il fait. $2$ $0.8$

Parcelle résiduelle

Enfin, nous pouvons décrire l'ajustement de plusieurs manières. L'ajustement est apparu dans les dernières lignes et colonnes de la matrice de nuage de points ( qv ) et mérite peut-être d'être examiné de plus près à ce point. En bas à gauche, les ajustements sont tracés sous forme de cercles et de flèches bleues ouvertes (représentant les résidus) les reliant aux données, représentées par des cercles rouges pleins. A droite, les sont représentés par des cercles noirs ouverts remplis de couleurs correspondant à leurs arguments; ceux-ci sont reliés par des flèches aux valeurs correspondantes de . Rappelons que chaque flèche représente une expansion par autour de l'origine, la rotation de $(w_j)$ $(z_j)$ $3/2$ $120$ degrés, et traduction par , plus cette erreur guassienne bivariée. $(-20, 5)$

Fit as transformation

Ces résultats, les graphiques et les graphiques de diagnostic suggèrent tous que la formule de régression complexe fonctionne correctement et réalise quelque chose de différent que les régressions linéaires séparées des parties réelle et imaginaire des variables.

Code

Le Rcode pour créer les données, les ajustements et les tracés apparaît ci-dessous. On notera que la solution réelle de est obtenu en une seule ligne de code. Un travail supplémentaire - mais pas trop - serait nécessaire pour obtenir la sortie habituelle des moindres carrés: la matrice variance-covariance de l'ajustement, les erreurs standard, les valeurs de p, etc. $\hat\beta$

#
# Synthesize data.
# (1) the independent variable `w`.
#
w.max <- 5 # Max extent of the independent values
w <- expand.grid(seq(-w.max,w.max), seq(-w.max,w.max))
w <- complex(real=w[[1]], imaginary=w[[2]])
w <- w[Mod(w) <= w.max]
n <- length(w)
#
# (2) the dependent variable `z`.
#
beta <- c(-20+5i, complex(argument=2*pi/3, modulus=3/2))
sigma <- 2; rho <- 0.8 # Parameters of the error distribution
library(MASS) #mvrnorm
set.seed(17)
e <- mvrnorm(n, c(0,0), matrix(c(1,rho,rho,1)*sigma^2, 2))
e <- complex(real=e[,1], imaginary=e[,2])
z <- as.vector((X <- cbind(rep(1,n), w)) %*% beta + e)
#
# Fit the models.
#
print(beta, digits=3)
print(beta.hat <- solve(Conj(t(X)) %*% X, Conj(t(X)) %*% z), digits=3)
print(beta.r <- coef(lm(Re(z) ~ Re(w) + Im(w))), digits=3)
print(beta.i <- coef(lm(Im(z) ~ Re(w) + Im(w))), digits=3)
#
# Show some diagnostics.
#
par(mfrow=c(1,2))
res <- as.vector(z - X %*% beta.hat)
fit <- z - res
s <- sqrt(Re(mean(Conj(res)*res)))
col <- hsv((Arg(res)/pi + 1)/2, .8, .9)
size <- Mod(res) / s
plot(res, pch=16, cex=size, col=col, main="Residuals")
plot(Re(fit), Im(fit), pch=16, cex = size, col=col,
     main="Residuals vs. Fitted")

plot(Re(c(z, fit)), Im(c(z, fit)), type="n",
     main="Residuals as Fit --> Data", xlab="Real", ylab="Imaginary")
points(Re(fit), Im(fit), col="Blue")
points(Re(z), Im(z), pch=16, col="Red")
arrows(Re(fit), Im(fit), Re(z), Im(z), col="Gray", length=0.1)

col.w <-  hsv((Arg(w)/pi + 1)/2, .8, .9)
plot(Re(c(w, z)), Im(c(w, z)), type="n",
     main="Fit as a Transformation", xlab="Real", ylab="Imaginary")
points(Re(w), Im(w), pch=16, col=col.w)
points(Re(w), Im(w))
points(Re(z), Im(z), pch=16, col=col.w)
arrows(Re(w), Im(w), Re(z), Im(z), col="#00000030", length=0.1)
#
# Display the data.
#
par(mfrow=c(1,1))
pairs(cbind(w.Re=Re(w), w.Im=Im(w), z.Re=Re(z), z.Im=Im(z),
            fit.Re=Re(fit), fit.Im=Im(fit)), cex=1/2)

whuber
la source

\hat{β}

$\hat{\beta}$

y

$y$

Si tout a été calculé correctement, la covariance sera toujours positive-définie. En particulier, cela implique que lorsque vous l'utilisez pour calculer la covariance de la partie réelle ou de la partie imaginaire d'une variable, vous obtiendrez un nombre positif, de sorte que tous les CI seront bien définis.

whuber

\hat{β}

$\hat{\beta}$

De plus, si lorsque je calcule des valeurs pour des statistiques de test, j'obtiens des nombres comme disons, 3 + .1 * i. Pour cela, je m'attendais à ce que le nombre n'ait pas de partie imaginaire. Est-ce normal? Ou un signe que je fais quelque chose de mal?

bill_e

Lorsque vous calculez des statistiques de test avec des nombres complexes, vous devez vous attendre à obtenir des résultats complexes! Si vous avez une raison mathématique pour laquelle la statistique doit être réelle, alors le calcul doit être erroné. Lorsque la partie imaginaire est vraiment minuscule par rapport à la partie réelle, il s'agit probablement d'une erreur en virgule flottante accumulée et il est généralement sûr de la tuer ( zapsmallin R). Sinon, c'est un signe que quelque chose cloche fondamentalement.

whuber

5

Après une longue session Google, j'ai trouvé des informations pertinentes pour comprendre le problème d'une manière alternative. Il s'avère que des problèmes similaires sont quelque peu courants dans le traitement statistique du signal. Au lieu de commencer par une vraisemblance gaussienne qui correspond à un moindre carré linéaire pour des données réelles, on commence par un:

http://en.wikipedia.org/wiki/Complex_normal_distribution

Cette page wikipedia donne un aperçu satisfaisant de cet objet.

$\hat{\beta}$

Une autre source que j'ai trouvée qui arrive à la même conclusion que whuber, mais explore d'autres estimateurs comme le maximum de vraisemblance est: «Estimations of the Improper Linear Regression Models», de Yan et al.

bill_e
la source

1

Alors que @whuber a une réponse magnifiquement illustrée et bien expliquée, je pense que c'est un modèle simplifié qui manque une partie de la puissance de l'espace complexe.

$w$ $\beta$ $x$

z = β_{0} + β_{1} w + ϵ

$z = \beta_0 + \beta_1 w + \epsilon$

$\epsilon$

Je suggère que la régression linéaire complexe soit définie comme suit:

z = β_{0} + β_{1} w + β_{2} \bar{w} + ϵ

$z = \beta_0 + \beta_1 w + \beta_2 \overline w + \epsilon$

Il existe deux différences majeures.

$\beta_2$

$\epsilon$

Pour en revenir au modèle réel, la solution des moindres carrés ordinaires sort minimisant la perte, qui est la log-vraisemblance négative. Pour une distribution normale, voici la parabole:

y = a x^{2} + c x + d .

$y = ax^2 + cx + d.$

$x = z - (\beta_0 + \beta_1 w)$ $a$ $c$ $d$

\begin{aligned} y = a {| x |}^{2} + ℜ (b x^{2} + c x) + d . \end{aligned}

$\begin{align} y = a{|x|}^2 + \Re({bx^2 + cx}) + d. \end{align}$

$c$ $d$ $a$ $b$ $b$ $\Re$

\begin{aligned} {[\begin{matrix} x - μ \\ \bar{x - μ} \end{matrix}]}^{H} {[\begin{matrix} s & u \\ \bar{u} & \bar{s} \end{matrix}]}^{- 1} [\begin{matrix} x - μ \\ \bar{x - μ} \end{matrix}] + d \end{aligned}

$\begin{align} {\begin{bmatrix}x-\mu \\ \overline{x-\mu}\end{bmatrix}}^H \begin{bmatrix}s & u \\ \overline{u} & \overline{s}\end{bmatrix}^{-1}\! \begin{bmatrix}x-\mu \\ \overline{x-\mu}\end{bmatrix} + d \end{align}$

s, u, μ, d

$s, u, \mu, d$

s

$s$

u

$u$

μ

$\mu$

Voici une image de la densité d'une distribution normale complexe:

$b$

Cela complique la régression bien que je sois presque sûr que la solution est encore analytique. Je l'ai résolu pour le cas d'une entrée, et je suis heureux de transcrire ma solution ici, mais j'ai le sentiment que whuber pourrait résoudre le cas général.

Neil G
la source

Merci pour cette contribution. Je ne le suis pas, cependant, car je ne sais pas (a) pourquoi vous introduisez un polynôme quadratique, (b) ce que vous entendez réellement par polynôme "correspondant", ou (c) quel modèle statistique vous ajustez. Pourriez-vous nous en dire plus?

whuber

@whuber, je l'ai réécrit comme modèle statistique. Veuillez me faire savoir si cela a du sens pour vous.

Neil G

z

$z$

w

$w$

\bar{w}

$\bar w$

ϵ

$\epsilon$

\Beta_{2}

$\Beta_2$

| x |^{2}

$|x|^2$

x^{2}

$x^2$

1

Ce problème est revenu sur Mathematica StackExchange et ma réponse / commentaire étendu est que l'excellente réponse de @whuber devrait être suivie.

Ma réponse ici est une tentative d'étendre un peu la réponse de @whuber en rendant la structure d'erreur un peu plus explicite. L'estimateur des moindres carrés proposé est celui que l'on utiliserait si la distribution d'erreur bivariée a une corrélation nulle entre les composantes réelle et imaginaire. (Mais les données générées ont une corrélation d'erreur de 0,8.)

$\rho=0$ $\rho\neq0$

$\rho=0$

$\rho$

$\gamma_0$ $\delta_0$ $\rho$ $\gamma_1$

Mon point dans tout cela est que le modèle qui convient doit être rendu complètement explicite et que les programmes d'algèbre symbolique peuvent aider à atténuer le désordre. (Et, bien sûr, les estimateurs du maximum de vraisemblance supposent une distribution normale bivariée que les estimateurs des moindres carrés ne supposent pas.)

Annexe: le code Mathematica complet

(* Predictor variable *)
w = {0 - 5 I, -3 - 4 I, -2 - 4 I, -1 - 4 I, 0 - 4 I, 1 - 4 I, 2 - 4 I,
    3 - 4 I, -4 - 3 I, -3 - 3 I, -2 - 3 I, -1 - 3 I, 0 - 3 I, 1 - 3 I,
    2 - 3 I, 3 - 3 I, 4 - 3 I, -4 - 2 I, -3 - 2 I, -2 - 2 I, -1 - 2 I,
    0 - 2 I, 1 - 2 I, 2 - 2 I, 3 - 2 I, 
   4 - 2 I, -4 - 1 I, -3 - 1 I, -2 - 1 I, -1 - 1 I, 0 - 1 I, 1 - 1 I, 
   2 - 1 I, 3 - 1 I, 
   4 - 1 I, -5 + 0 I, -4 + 0 I, -3 + 0 I, -2 + 0 I, -1 + 0 I, 0 + 0 I,
    1 + 0 I, 2 + 0 I, 3 + 0 I, 4 + 0 I, 
   5 + 0 I, -4 + 1 I, -3 + 1 I, -2 + 1 I, -1 + 1 I, 0 + 1 I, 1 + 1 I, 
   2 + 1 I, 3 + 1 I, 4 + 1 I, -4 + 2 I, -3 + 2 I, -2 + 2 I, -1 + 2 I, 
   0 + 2 I, 1 + 2 I, 2 + 2 I, 3 + 2 I, 
   4 + 2 I, -4 + 3 I, -3 + 3 I, -2 + 3 I, -1 + 3 I, 0 + 3 I, 1 + 3 I, 
   2 + 3 I, 3 + 3 I, 4 + 3 I, -3 + 4 I, -2 + 4 I, -1 + 4 I, 0 + 4 I, 
   1 + 4 I, 2 + 4 I, 3 + 4 I, 0 + 5 I};
(* Add in a "1" for the intercept *)
w1 = Transpose[{ConstantArray[1 + 0 I, Length[w]], w}];

z = {-15.83651 + 7.23001 I, -13.45474 + 4.70158 I, -13.63353 + 
    4.84748 I, -14.79109 + 4.33689 I, -13.63202 + 
    9.75805 I, -16.42506 + 9.54179 I, -14.54613 + 
    12.53215 I, -13.55975 + 14.91680 I, -12.64551 + 
    2.56503 I, -13.55825 + 4.44933 I, -11.28259 + 
    5.81240 I, -14.14497 + 7.18378 I, -13.45621 + 
    9.51873 I, -16.21694 + 8.62619 I, -14.95755 + 
    13.24094 I, -17.74017 + 10.32501 I, -17.23451 + 
    13.75955 I, -14.31768 + 1.82437 I, -13.68003 + 
    3.50632 I, -14.72750 + 5.13178 I, -15.00054 + 
    6.13389 I, -19.85013 + 6.36008 I, -19.79806 + 
    6.70061 I, -14.87031 + 11.41705 I, -21.51244 + 
    9.99690 I, -18.78360 + 14.47913 I, -15.19441 + 
    0.49289 I, -17.26867 + 3.65427 I, -16.34927 + 
    3.75119 I, -18.58678 + 2.38690 I, -20.11586 + 
    2.69634 I, -22.05726 + 6.01176 I, -22.94071 + 
    7.75243 I, -28.01594 + 3.21750 I, -24.60006 + 
    8.46907 I, -16.78006 - 2.66809 I, -18.23789 - 
    1.90286 I, -20.28243 + 0.47875 I, -18.37027 + 
    2.46888 I, -21.29372 + 3.40504 I, -19.80125 + 
    5.76661 I, -21.28269 + 5.57369 I, -22.05546 + 
    7.37060 I, -18.92492 + 10.18391 I, -18.13950 + 
    12.51550 I, -22.34471 + 10.37145 I, -15.05198 + 
    2.45401 I, -19.34279 - 0.23179 I, -17.37708 + 
    1.29222 I, -21.34378 - 0.00729 I, -20.84346 + 
    4.99178 I, -18.01642 + 10.78440 I, -23.08955 + 
    9.22452 I, -23.21163 + 7.69873 I, -26.54236 + 
    8.53687 I, -16.19653 - 0.36781 I, -23.49027 - 
    2.47554 I, -21.39397 - 0.05865 I, -20.02732 + 
    4.10250 I, -18.14814 + 7.36346 I, -23.70820 + 
    5.27508 I, -25.31022 + 4.32939 I, -24.04835 + 
    7.83235 I, -26.43708 + 6.19259 I, -21.58159 - 
    0.96734 I, -21.15339 - 1.06770 I, -21.88608 - 
    1.66252 I, -22.26280 + 4.00421 I, -22.37417 + 
    4.71425 I, -27.54631 + 4.83841 I, -24.39734 + 
    6.47424 I, -30.37850 + 4.07676 I, -30.30331 + 
    5.41201 I, -28.99194 - 8.45105 I, -24.05801 + 
    0.35091 I, -24.43580 - 0.69305 I, -29.71399 - 
    2.71735 I, -26.30489 + 4.93457 I, -27.16450 + 
    2.63608 I, -23.40265 + 8.76427 I, -29.56214 - 2.69087 I};

(* whuber 's least squares estimates *)
{a, b} = Inverse[ConjugateTranspose[w1].w1].ConjugateTranspose[w1].z
(* {-20.0172+5.00968 \[ImaginaryI],-0.830797+1.37827 \[ImaginaryI]} *)

(* Break up into the real and imaginary components *)
x = Re[z];
y = Im[z];
u = Re[w];
v = Im[w];
n = Length[z]; (* Sample size *)

(* Construct the real and imaginary components of the model *)
(* This is the messy part you probably don't want to do too often with paper and pencil *)
model = \[Gamma]0 + I \[Delta]0 + (\[Gamma]1 + I \[Delta]1) (u + I v);
modelR = Table[
   Re[ComplexExpand[model[[j]]]] /. Im[h_] -> 0 /. Re[h_] -> h, {j, n}];
(* \[Gamma]0+u \[Gamma]1-v \[Delta]1 *)
modelI = Table[
   Im[ComplexExpand[model[[j]]]] /. Im[h_] -> 0 /. Re[h_] -> h, {j, n}];
(* v \[Gamma]1+\[Delta]0+u \[Delta]1 *)

(* Construct the log of the likelihood as we are estimating the parameters associated with a bivariate normal distribution *)
logL = LogLikelihood[
   BinormalDistribution[{0, 0}, {\[Sigma]1, \[Sigma]2}, \[Rho]],
   Transpose[{x - modelR, y - modelI}]];

mle0 = FindMaximum[{logL /. {\[Rho] -> 
      0, \[Sigma]1 -> \[Sigma], \[Sigma]2 -> \[Sigma]}, \[Sigma] > 
    0}, {\[Gamma]0, \[Delta]0, \[Gamma]1, \[Delta]1, \[Sigma]}]
(* {-357.626,{\[Gamma]0\[Rule]-20.0172,\[Delta]0\[Rule]5.00968,\[Gamma]1\[Rule]-0.830797,\[Delta]1\[Rule]1.37827,\[Sigma]\[Rule]2.20038}} *)

(* Now suppose we don't want to restrict \[Rho]=0 *)
mle1 = FindMaximum[{logL /. {\[Sigma]1 -> \[Sigma], \[Sigma]2 -> \[Sigma]}, \[Sigma] > 0 && -1 < \[Rho] < 
     1}, {\[Gamma]0, \[Delta]0, \[Gamma]1, \[Delta]1, \[Sigma], \[Rho]}]
(* {-315.313,{\[Gamma]0\[Rule]-20.0172,\[Delta]0\[Rule]5.00968,\[Gamma]1\[Rule]-0.763237,\[Delta]1\[Rule]1.30859,\[Sigma]\[Rule]2.21424,\[Rho]\[Rule]0.810525}} *)

JimB
la source

Analyse avec des données complexes, quelque chose de différent?

Réponses:

Sommaire

Interprétation

Les moindres carrés ordinaires

Exemple

Code