{VERSION 4 0 "IBM INTEL NT" "4.0" } {USTYLETAB {CSTYLE "Maple Input" -1 0 "Courier" 0 1 255 0 0 1 0 1 0 0 1 0 0 0 0 1 }{CSTYLE "2D Math" -1 2 "Times" 0 1 0 0 0 0 0 0 2 0 0 0 0 0 0 1 }{CSTYLE "2D Comment" 2 18 "" 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 1 } {CSTYLE "" -1 256 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 257 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 258 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 259 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 260 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 } {CSTYLE "" -1 261 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 262 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 263 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 264 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 265 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 } {CSTYLE "" -1 266 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 267 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 268 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 269 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 270 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 } {CSTYLE "" -1 271 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 272 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 273 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 274 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 275 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 } {CSTYLE "" -1 276 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 277 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 278 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 279 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 280 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 } {CSTYLE "" -1 281 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 282 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 283 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 284 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 285 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 } {CSTYLE "" -1 286 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 287 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 288 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 289 "" 1 12 0 0 0 0 0 2 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 290 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 } {CSTYLE "" -1 291 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 292 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 293 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 294 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 295 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 } {CSTYLE "" -1 296 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 297 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 298 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 299 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 300 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 } {CSTYLE "" -1 301 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 302 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 303 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 304 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 305 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 } {CSTYLE "" -1 306 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 307 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 308 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 309 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 310 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 } {CSTYLE "" -1 311 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 312 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 313 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 314 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 315 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 } {CSTYLE "" -1 316 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 317 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 318 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 319 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 320 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 } {CSTYLE "" -1 321 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 322 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 323 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 324 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 325 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 } {CSTYLE "" -1 326 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 327 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 328 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 329 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 330 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 } {CSTYLE "" -1 331 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 332 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 333 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 334 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 335 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 } {CSTYLE "" -1 336 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 337 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 338 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 339 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 340 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 } {CSTYLE "" -1 341 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 342 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 343 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 344 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 345 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 } {CSTYLE "" -1 346 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 347 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 348 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 349 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 350 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 } {CSTYLE "" -1 351 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 352 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 353 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 354 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 355 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 } {CSTYLE "" -1 356 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 357 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 358 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 359 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 360 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 } {CSTYLE "" -1 361 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 362 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 363 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 364 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 365 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 } {CSTYLE "" -1 366 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 367 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 368 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 369 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 370 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 } {CSTYLE "" -1 371 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 372 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 373 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 374 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 375 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 } {CSTYLE "" -1 376 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 377 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 378 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 379 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 380 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 } {CSTYLE "" -1 381 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 382 "" 1 24 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 383 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 384 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 385 "" 1 24 0 0 0 0 0 1 0 0 0 0 0 0 0 1 } {CSTYLE "" -1 386 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 387 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 388 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 389 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 390 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 } {CSTYLE "" -1 391 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 392 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 393 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 394 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 395 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 } {CSTYLE "" -1 396 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 397 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 398 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 399 "Verdana" 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 400 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 401 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 402 "Verdana" 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 403 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 404 "Verdana" 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 405 "Verdana" 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 406 "Verdana" 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 407 "Verdana" 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 1 } {CSTYLE "" -1 408 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 409 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 410 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 411 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 412 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 } {CSTYLE "" -1 413 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 414 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 415 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 416 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 417 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 } {CSTYLE "" -1 418 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 419 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 420 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 421 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 422 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 } {CSTYLE "" -1 423 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 424 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 425 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 426 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 427 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 } {CSTYLE "" -1 428 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 429 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 430 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 431 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 432 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 } {CSTYLE "" -1 433 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 434 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 435 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 436 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 437 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 } {CSTYLE "" -1 438 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 439 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 440 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 441 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 442 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 } {CSTYLE "" -1 443 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 444 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 445 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 446 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 447 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 } {CSTYLE "" -1 448 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 449 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 450 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 451 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 452 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 } {CSTYLE "" -1 453 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 454 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 455 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 456 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 457 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 } {CSTYLE "" -1 458 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 459 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 460 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 461 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 462 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 } {CSTYLE "" -1 463 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 464 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 465 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 466 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 467 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 } {CSTYLE "" -1 468 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 469 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 470 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 471 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 472 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 } {CSTYLE "" -1 473 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 474 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 475 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 476 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 477 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 } {CSTYLE "" -1 478 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 479 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 480 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 481 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 482 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 } {CSTYLE "" -1 483 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 484 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 485 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 486 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 487 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 } {CSTYLE "" -1 488 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 489 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 490 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 491 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 492 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 } {CSTYLE "" -1 493 "Verdana" 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 1 }{CSTYLE " " -1 494 "Verdana" 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 495 "Verdana" 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 496 "Verdana " 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 497 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 498 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 499 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 500 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 501 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 502 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 503 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 } {CSTYLE "" -1 504 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 505 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 506 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 507 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 508 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 } {CSTYLE "" -1 509 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 510 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 511 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 512 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 513 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 } {CSTYLE "" -1 514 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{PSTYLE "Normal " -1 0 1 {CSTYLE "" -1 -1 "" 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 1 }0 0 0 -1 -1 -1 0 0 0 0 0 0 -1 0 }{PSTYLE "Heading 1" 0 3 1 {CSTYLE "" -1 -1 "" 1 18 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }1 0 0 0 6 6 0 0 0 0 0 0 -1 0 } {PSTYLE "R3 Font 0" -1 256 1 {CSTYLE "" -1 -1 "Times" 0 14 0 0 0 1 2 1 2 0 0 0 0 0 0 1 }0 0 0 -1 -1 -1 0 0 0 0 0 0 -1 0 }{PSTYLE "R3 Font 2 " -1 257 1 {CSTYLE "" -1 -1 "Courier" 0 11 0 0 0 1 2 1 2 0 0 0 0 0 0 1 }0 0 0 -1 -1 -1 0 0 0 0 0 0 -1 0 }{PSTYLE "" 0 258 1 {CSTYLE "" -1 -1 "" 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 1 }3 0 0 -1 -1 -1 0 0 0 0 0 0 -1 0 }{PSTYLE "" 0 259 1 {CSTYLE "" -1 -1 "" 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 1 }3 0 0 -1 -1 -1 0 0 0 0 0 0 -1 0 }{PSTYLE "" 0 260 1 {CSTYLE "" -1 -1 "" 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 1 }3 0 0 -1 -1 -1 0 0 0 0 0 0 -1 0 }{PSTYLE "" 0 261 1 {CSTYLE "" -1 -1 "" 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 1 }3 0 0 -1 -1 -1 0 0 0 0 0 0 -1 0 }{PSTYLE "" 0 262 1 {CSTYLE " " -1 -1 "" 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 1 }3 0 0 -1 -1 -1 0 0 0 0 0 0 -1 0 }{PSTYLE "" 0 263 1 {CSTYLE "" -1 -1 "" 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 1 }3 0 0 -1 -1 -1 0 0 0 0 0 0 -1 0 }{PSTYLE "" 0 264 1 {CSTYLE "" -1 -1 "" 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 1 }3 0 0 -1 -1 -1 0 0 0 0 0 0 -1 0 }{PSTYLE "" 0 265 1 {CSTYLE "" -1 -1 "" 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 1 }3 0 0 -1 -1 -1 0 0 0 0 0 0 -1 0 }{PSTYLE "" 0 266 1 {CSTYLE "" -1 -1 "" 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 1 }3 0 0 -1 -1 -1 0 0 0 0 0 0 -1 0 }{PSTYLE "" 0 267 1 {CSTYLE "" -1 -1 "" 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 1 }3 0 0 -1 -1 -1 0 0 0 0 0 0 -1 0 }{PSTYLE "" 0 268 1 {CSTYLE "" -1 -1 "" 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 1 }3 0 0 -1 -1 -1 0 0 0 0 0 0 -1 0 }{PSTYLE "" 0 269 1 {CSTYLE "" -1 -1 "" 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 1 }3 0 0 -1 -1 -1 0 0 0 0 0 0 -1 0 }{PSTYLE " " 0 270 1 {CSTYLE "" -1 -1 "" 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 1 }3 0 0 -1 -1 -1 0 0 0 0 0 0 -1 0 }{PSTYLE "" 0 271 1 {CSTYLE "" -1 -1 "" 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 1 }3 0 0 -1 -1 -1 0 0 0 0 0 0 -1 0 }{PSTYLE "" 0 272 1 {CSTYLE "" -1 -1 "" 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 1 }3 0 0 -1 -1 -1 0 0 0 0 0 0 -1 0 }{PSTYLE "" 0 273 1 {CSTYLE "" -1 -1 "" 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 1 }3 0 0 -1 -1 -1 0 0 0 0 0 0 -1 0 }{PSTYLE "" 0 274 1 {CSTYLE "" -1 -1 "" 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 1 }3 0 0 -1 -1 -1 0 0 0 0 0 0 -1 0 }{PSTYLE "" 0 275 1 {CSTYLE "" -1 -1 "" 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 1 }3 0 0 -1 -1 -1 0 0 0 0 0 0 -1 0 }} {SECT 0 {EXCHG {PARA 0 "" 0 "" {TEXT -1 22 "Unit 6: Matrix Algebra" }} {PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 0 "" 0 "" {TEXT -1 33 "Chapter 3 4: Matrix Factorizations" }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 0 " " 0 "" {TEXT -1 76 "Section 34.6: minimum-length least-squares soluti ons, and the pseudoinverse" }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}}{SECT 1 {PARA 3 "" 0 "" {TEXT -1 9 "Copyright" }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{EXCHG {PARA 3 "" 0 "" {TEXT 289 48 "Copyright * 2001 by Addison \+ Wesley Longman, Inc." }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 0 "" 0 " " {TEXT -1 302 "All rights reserved. No part of this publication may \+ be reproduced, stored in a retrieval system, or transmitted, in any fo rm or by any means, electronic, mechanical, photocopying, recording, o r otherwise, without the prior written permission of the publisher. Pr inted in the United States of America." }}{PARA 0 "" 0 "" {TEXT -1 0 " " }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" }}}}{SECT 1 {PARA 3 " " 0 "" {TEXT -1 15 "Initializations" }}{PARA 0 "" 0 "" {TEXT -1 0 "" } }{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 8 "restart;" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 61 "with(linalg):\nwith(plots):\nwith(plottools ):\nread(`pvac.txt`):" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 29 "al ias(c[1]=_t[1],c[2]=_t[2]);" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 17 "f := rand(-9..9):" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" } }}}{SECT 1 {PARA 3 "" 0 "" {TEXT -1 34 "Least Squares Problem: Formula tion" }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{EXCHG {PARA 0 "" 0 "" {TEXT -1 307 "In this section we again solve a least squares problem as an e xample of an over-determined system. However, we will look at the gen eral case where the rank of the coefficient matrix is less than maxima l. Since the system is over-determined, there are more rows than colu mns in the system matrix, that is, " }{XPPEDIT 18 0 "r*`>`*c" "6#*(%\" rG\"\"\"%\">GF%%\"cGF%" }{TEXT -1 39 ". If the rank of the system mat rix is " }{XPPEDIT 18 0 "k " 0 "" {MPLTEXT 1 0 101 "_seed := 125:\na1 := randmatrix(5,2,entries=f):\na2 := randma trix(2,4,entries=f):\nA := evalm(a1 &* a2);" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 23 "To generat e the vector " }{TEXT 312 1 "y" }{TEXT -1 39 ", execute the following \+ Maple commands." }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}}{EXCHG {PARA 0 "> \+ " 0 "" {MPLTEXT 1 0 52 "_seed := 718377433241:\ny := randvector(5,entr ies=f);" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 11 "The matrix " }{TEXT 345 1 "A" }{TEXT -1 46 " indeed has rank 2, as we saw in Section 34.5." }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 8 "rank(A);" }} }{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" }}}}{SECT 1 {PARA 3 "" 0 " " {TEXT -1 31 "Least Squares Problem: Solution" }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{EXCHG {PARA 0 "" 0 "" {TEXT -1 23 "The third row of r ref([" }{TEXT 349 1 "A" }{TEXT -1 2 ", " }{TEXT 350 1 "y" }{TEXT -1 58 "]) given in Table 34.4 indicates an inconsistent equation " } {XPPEDIT 18 0 "0=1" "6#/\"\"!\"\"\"" }{TEXT -1 74 ", so we know there \+ is no solution to the over-determined system of linear " }{TEXT 346 1 "A" }{TEXT 347 1 "x" }{TEXT -1 3 " = " }{TEXT 348 1 "y" }{TEXT -1 54 " . Hence, we seek a least squares \"best fit\" solution." }}{PARA 0 " " 0 "" {TEXT -1 0 "" }}{PARA 0 "" 0 "" {TEXT -1 34 "The normal equatio ns are given by " }{XPPEDIT 18 0 "A^T*A" "6#*&)%\"AG%\"TG\"\"\"F%F'" } {TEXT 351 1 "x" }{TEXT -1 3 " = " }{XPPEDIT 18 0 "A^T" "6#)%\"AG%\"TG " }{TEXT 352 1 "y" }{TEXT -1 8 ", where " }{XPPEDIT 18 0 "A^T*A" "6#*& )%\"AG%\"TG\"\"\"F%F'" }{TEXT -1 2 ", " }{TEXT 353 1 "b" }{TEXT -1 3 " = " }{XPPEDIT 18 0 "A^T" "6#)%\"AG%\"TG" }{TEXT 354 1 "y" }{TEXT -1 19 ", and the solution " }{TEXT 355 1 "x" }{TEXT -1 3 " = " }{TEXT 356 1 "X" }{TEXT -1 5 ", are" }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 262 "" 0 "" {TEXT -1 1 " " }{XPPEDIT 18 0 "A^T*A=MATRIX([[13286, -1377 6, 1162, -938], [-13776, 22432, -5464, 6528], [1162, -5464, 2328, -298 6], [-938, 6528, -2986, 3854]])" "6#/*&)%\"AG%\"TG\"\"\"F&F(-%'MATRIXG 6#7&7&\"&'G8,$\"&wP\"!\"\"\"%i6,$\"$Q*F17&,$F0F1\"&KC#,$\"%kaF1\"%Gl7& F2,$F9F1\"%GB,$\"%')HF17&,$F4F1F:,$F?F1\"%aQ" }{TEXT -1 5 ", " } {TEXT 357 1 "b" }{TEXT -1 3 " = " }{XPPEDIT 18 0 "MATRIX([[-254], [328 ], [-56], [62]])" "6#-%'MATRIXG6#7&7#,$\"$a#!\"\"7#\"$G$7#,$\"#cF*7#\" #i" }{TEXT -1 5 ", " }{TEXT 358 1 "X" }{TEXT -1 3 " = " }{XPPEDIT 18 0 "MATRIX([[c[1]], [c[2]], [-30*c[1]+28*c[2]-7673/13979], [-23*c[1] +20*c[2]-5720/13979]])" "6#-%'MATRIXG6#7&7#&%\"cG6#\"\"\"7#&F)6#\"\"#7 #,(*&\"#IF+&F)6#F+F+!\"\"*&\"#GF+&F)6#F/F+F+*&\"%twF+\"&zR\"F6F67#,(*& \"#BF+&F)6#F+F+F6*&\"#?F+&F)6#F/F+F+*&\"%?dF+F=F6F6" }{TEXT -1 1 " " } }{PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 0 "" 0 "" {TEXT -1 27 "The leas t-squares solution " }{TEXT 359 1 "X" }{TEXT -1 38 " contains two free parameters because " }{XPPEDIT 18 0 "A^T*A" "6#*&)%\"AG%\"TG\"\"\"F%F '" }{TEXT -1 7 ", like " }{TEXT 360 1 "A" }{TEXT -1 230 ", has rank 2. There are only two independent equations in the set of four normal e quations. Note, however, that the four normal equations are now consi stent, but under-determined. There is a solution, but not a unique so lution." }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 0 "" 0 "" {TEXT -1 28 "We now write the solution as" }}{PARA 263 "" 0 "" {TEXT -1 1 " " } {TEXT 361 1 "X" }{TEXT -1 3 " = " }{TEXT 362 1 "v" }{XPPEDIT 18 0 "``[ 1]" "6#&%!G6#\"\"\"" }{TEXT -1 3 " + " }{XPPEDIT 18 0 "c[1]" "6#&%\"cG 6#\"\"\"" }{TEXT 363 1 "v" }{XPPEDIT 18 0 "``[2]" "6#&%!G6#\"\"#" } {TEXT -1 3 " + " }{XPPEDIT 18 0 "c[2]" "6#&%\"cG6#\"\"#" }{TEXT 364 1 "v" }{XPPEDIT 18 0 "``[3]" "6#&%!G6#\"\"$" }{TEXT -1 1 " " }}{PARA 0 " " 0 "" {TEXT -1 6 "where " }{TEXT 365 1 "v" }{XPPEDIT 18 0 "``[k],k=1, 2,3" "6&&%!G6#%\"kG/F&\"\"\"\"\"#\"\"$" }{TEXT -1 14 ", are given by" }}{PARA 264 "" 0 "" {TEXT -1 1 " " }{TEXT 366 1 "v" }{XPPEDIT 18 0 "`` [1]" "6#&%!G6#\"\"\"" }{TEXT -1 3 " = " }{XPPEDIT 18 0 "MATRIX([[0], [ 0], [-7673/13979], [-5720/13979]])" "6#-%'MATRIXG6#7&7#\"\"!7#F(7#,$*& \"%tw\"\"\"\"&zR\"!\"\"F07#,$*&\"%?dF.F/F0F0" }{TEXT -1 5 ", " } {TEXT 367 1 "v" }{XPPEDIT 18 0 "``[2]" "6#&%!G6#\"\"#" }{TEXT -1 3 " = " }{XPPEDIT 18 0 "MATRIX([[1], [0], [-30], [-23]])" "6#-%'MATRIXG6#7& 7#\"\"\"7#\"\"!7#,$\"#I!\"\"7#,$\"#BF." }{TEXT -1 5 ", " }{TEXT 368 1 "v" }{XPPEDIT 18 0 "``[3]" "6#&%!G6#\"\"$" }{TEXT -1 3 " = " } {XPPEDIT 18 0 "MATRIX([[0], [1], [28], [20]])" "6#-%'MATRIXG6#7&7#\"\" !7#\"\"\"7#\"#G7#\"#?" }{TEXT -1 1 " " }}{PARA 0 "" 0 "" {TEXT -1 0 " " }}{PARA 0 "" 0 "" {TEXT -1 12 "The vectors " }{TEXT 369 1 "v" } {XPPEDIT 18 0 "``[2]" "6#&%!G6#\"\"#" }{TEXT -1 5 " and " }{TEXT 370 1 "v" }{XPPEDIT 18 0 "``[3]" "6#&%!G6#\"\"$" }{TEXT -1 26 " are in the null space of " }{TEXT 372 1 "A" }{TEXT -1 9 " because " }{TEXT 371 1 "A" }{TEXT 373 1 "v" }{XPPEDIT 18 0 "``[k]" "6#&%!G6#%\"kG" }{TEXT -1 3 " = " }{TEXT 379 1 "0" }{TEXT -1 2 ", " }{XPPEDIT 18 0 "k=2,3" "6 $/%\"kG\"\"#\"\"$" }{TEXT -1 3 ". " }}{PARA 0 "" 0 "" {TEXT -1 0 "" } }{PARA 0 "" 0 "" {TEXT -1 11 "The vector " }{TEXT 375 1 "v" }{XPPEDIT 18 0 "``[1]" "6#&%!G6#\"\"\"" }{TEXT -1 4 " is " }{TEXT 374 3 "not" } {TEXT -1 15 " orthogonal to " }{TEXT 376 1 "v" }{XPPEDIT 18 0 "``[2]" "6#&%!G6#\"\"#" }{TEXT -1 5 " and " }{TEXT 377 1 "v" }{XPPEDIT 18 0 "` `[3]" "6#&%!G6#\"\"$" }{TEXT -1 8 " because" }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 265 "" 0 "" {TEXT -1 1 " " }{TEXT 380 1 "v" }{XPPEDIT 18 0 "``[1]" "6#&%!G6#\"\"\"" }{TEXT -1 1 " " }{TEXT 382 1 "." }{TEXT -1 1 " " }{TEXT 381 1 "v" }{XPPEDIT 18 0 "``[2]" "6#&%!G6#\"\"#" } {TEXT -1 3 " = " }{XPPEDIT 18 0 "361750/13979<>0" "6#0*&\"']0" "6#0,$*&\"'W#H$\"\"\"\"&zR\"!\" \"F)\"\"!" }{TEXT -1 1 " " }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 0 " " 0 "" {TEXT -1 6 "Thus, " }{TEXT 386 1 "v" }{XPPEDIT 18 0 "``[1]" "6# &%!G6#\"\"\"" }{TEXT -1 41 " is not orthogonal to the null space, so \+ " }{TEXT 387 1 "v" }{XPPEDIT 18 0 "``[1]" "6#&%!G6#\"\"\"" }{TEXT -1 28 " is not in the row space of " }{TEXT 388 1 "A" }{TEXT -1 14 ". Th e vector " }{TEXT 389 1 "v" }{XPPEDIT 18 0 "``[1]" "6#&%!G6#\"\"\"" } {TEXT -1 37 " has a component in the row space of " }{TEXT 390 1 "A" } {TEXT -1 51 ", but it also has a component in the null space of " } {TEXT 391 1 "A" }{TEXT -1 1 "." }}{PARA 0 "" 0 "" {TEXT -1 0 "" }} {PARA 0 "" 0 "" {TEXT -1 118 "Before we move on to the next conceptual segment of our discussion, we pause to implement these calculations i n Maple." }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 0 "" 0 "" {TEXT -1 61 "There are several ways to approach solving the linear system " } {TEXT 392 1 "A" }{TEXT 393 1 "x" }{TEXT -1 3 " = " }{TEXT 394 1 "y" } {TEXT -1 24 ". For example, Maple's " }{TEXT 256 8 "linsolve" }{TEXT -1 43 " command can be used. Here, the outcome is" }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 14 "linsolve(A, y);" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" }}}{EXCHG {PARA 0 " " 0 "" {TEXT -1 164 "The \"null\" response suggests there is no soluti on to the set of equations. We can verify that the equations are inco nsistent by row-reducing the augmented matrix [" }{TEXT 378 1 "A" } {TEXT -1 2 ", " }{TEXT 313 1 "y" }{TEXT -1 34 "] to its reduced row-ec helon form." }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 19 "rref(augment(A,y));" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 202 "The third ro w indicates an inconsistent equation 0 = 1, so we know there is no sol ution to this over-determined system of linear equations. That is why we will seek a least squares \"best fit\" solution." }}{PARA 0 "" 0 " " {TEXT -1 0 "" }}{PARA 0 "" 0 "" {TEXT -1 56 "First, we try creating \+ and solving the normal equations " }{TEXT 395 1 "A" }{XPPEDIT 18 0 "`` ^T" "6#)%!G%\"TG" }{TEXT 396 1 "A" }{TEXT -1 1 " " }{TEXT 257 1 "x" } {TEXT -1 3 " = " }{TEXT 397 1 "A" }{XPPEDIT 18 0 "``^T" "6#)%!G%\"TG" }{TEXT -1 1 " " }{TEXT 314 1 "y" }{TEXT -1 3 ". " }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 0 "" 0 "" {TEXT -1 6 "Since " }{TEXT 398 1 "A" } {TEXT -1 6 " is 5 " }{TEXT 399 1 "x" }{TEXT -1 4 " 4, " }{TEXT 400 1 " A" }{XPPEDIT 18 0 "``^T" "6#)%!G%\"TG" }{TEXT 401 1 "A" }{TEXT -1 10 " is the 4 " }{TEXT 402 1 "x" }{TEXT -1 10 " 4 matrix " }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 31 "ATA := e valm(transpose(A) &*A);" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 " " }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 15 "and the vector " }{TEXT 316 1 "b" }{TEXT -1 3 " = " }{TEXT 403 1 "A" }{XPPEDIT 18 0 "``^T" "6#)%!G %\"TG" }{TEXT 315 1 "y" }{TEXT -1 18 ", of dimension (4 " }{TEXT 404 1 "x" }{TEXT -1 6 " 5)(5 " }{TEXT 405 1 "x" }{TEXT -1 8 " 1) = 4 " } {TEXT 406 1 "x" }{TEXT -1 6 " 1, is" }}{PARA 0 "" 0 "" {TEXT -1 0 "" } }}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 30 "b := evalm(transpose(A) &* y);" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 61 "Hence, the normal equations form a \"square\" sys tem for which " }{TEXT 290 8 "linsolve" }{TEXT -1 20 " yields the solu tion" }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 21 "X := linsolve(ATA,b);" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 84 "The solution \+ of the normal equations contains two free parameters. The system is 4 " }{TEXT 407 1 "x" }{TEXT -1 13 " 4, and like " }{TEXT 408 1 "A" } {TEXT -1 14 ", the rank of " }{TEXT 409 1 "A" }{XPPEDIT 18 0 "``^T" "6 #)%!G%\"TG" }{TEXT 410 1 "A" }{TEXT -1 8 " is also" }}{PARA 0 "" 0 "" {TEXT -1 23 "of rank two, as seen by" }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 10 "rank(ATA);" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 259 "That is why there are two free parameters. There are only two independent equations in the set of four normal equations. Note, how ever, that the four normal equations are now consistent, but under-det ermined. There is a solution, but not a unique solution." }}{PARA 0 " " 0 "" {TEXT -1 0 "" }}{PARA 0 "" 0 "" {TEXT -1 28 "Next, we split the solution " }{TEXT 330 1 "X" }{TEXT -1 15 " into a vector " }{TEXT 270 1 "v" }{XPPEDIT 18 0 "``[1]" "6#&%!G6#\"\"\"" }{TEXT -1 14 ", and \+ vectors " }{TEXT 271 1 "v" }{XPPEDIT 18 0 "``[2]" "6#&%!G6#\"\"#" } {TEXT -1 5 " and " }{TEXT 272 1 "v" }{XPPEDIT 18 0 "``[3]" "6#&%!G6#\" \"$" }{TEXT -1 22 " in the null space of " }{TEXT 411 1 "A" }{TEXT -1 50 ". The splitting is accomplished by first setting " }{XPPEDIT 18 0 "c[1]=c[2]" "6#/&%\"cG6#\"\"\"&F%6#\"\"#" }{TEXT -1 78 " = 0, thereb y suppressing the part of the solution lying in the null space of " } {TEXT 412 1 "A" }{TEXT -1 42 ". Consequently, what's then left must b e " }{TEXT 291 1 "v" }{XPPEDIT 18 0 "``[1]" "6#&%!G6#\"\"\"" }{TEXT -1 1 "." }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 34 "v[1] := subs(c[1]=0,c[2]=0,op(X));" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 28 "T hen, after subtracting out " }{TEXT 273 1 "v" }{XPPEDIT 18 0 "``[1]" " 6#&%!G6#\"\"\"" }{TEXT -1 6 " from " }{TEXT 331 1 "X" }{TEXT -1 11 ", \+ obtaining" }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}}{EXCHG {PARA 0 "> " 0 " " {MPLTEXT 1 0 20 "vv := evalm(X-v[1]);" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 28 "we alternativ ely set either " }{XPPEDIT 18 0 "c[1]" "6#&%\"cG6#\"\"\"" }{TEXT -1 4 " or " }{XPPEDIT 18 0 "c[2]" "6#&%\"cG6#\"\"#" }{TEXT -1 38 " to zero \+ in the difference to extract " }{TEXT 292 1 "v" }{XPPEDIT 18 0 "``[2] " "6#&%!G6#\"\"#" }{TEXT -1 5 " and " }{TEXT 293 1 "v" }{XPPEDIT 18 0 "``[3]" "6#&%!G6#\"\"$" }{TEXT -1 52 ", the individual vectors spannin g the null space of " }{TEXT 413 1 "A" }{TEXT -1 1 "." }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 71 "v[2] := \+ subs(c[1]=1,c[2]=0,op(vv));\nv[3] := subs(c[1]=0,c[2]=1,op(vv));" }}} {EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 20 "To demonstrate that " }{TEXT 258 1 "v" }{XPPEDIT 18 0 "`` [2]" "6#&%!G6#\"\"#" }{TEXT -1 5 " and " }{TEXT 259 1 "v" }{XPPEDIT 18 0 "``[3]" "6#&%!G6#\"\"$" }{TEXT -1 26 " are in the null space of \+ " }{TEXT 414 1 "A" }{TEXT -1 7 ", show " }{TEXT 415 1 "A" }{TEXT 317 1 "v" }{XPPEDIT 18 0 "``[2]" "6#&%!G6#\"\"#" }{TEXT -1 3 " = " }{TEXT 260 1 "0" }{TEXT -1 5 " and " }{TEXT 416 1 "A" }{TEXT 261 1 "v" } {XPPEDIT 18 0 "``[3]" "6#&%!G6#\"\"$" }{TEXT -1 3 " = " }{TEXT 262 1 " 0" }{TEXT -1 10 ". Hence, " }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}} {EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 41 "print(evalm(A &* v[2]),evalm (A &* v[3]));" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" }}} {EXCHG {PARA 0 "" 0 "" {TEXT -1 81 "Next, show that the same least-squ ares solution is generated by applying Maple's " }{TEXT 263 9 "leastsq rs" }{TEXT -1 32 " command to the original system " }{TEXT 417 1 "A" } {TEXT 264 1 "x" }{TEXT -1 3 " = " }{TEXT 318 1 "y" }{TEXT -1 1 "." }} {PARA 0 "" 0 "" {TEXT -1 0 "" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 32 "Z := leastsqrs(A,y):\nprint(Z,X);" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 38 "Observe, if y ou will, that the vector " }{TEXT 274 1 "v" }{XPPEDIT 18 0 "``[1]" "6# &%!G6#\"\"\"" }{TEXT -1 4 " is " }{TEXT 287 3 "not" }{TEXT -1 15 " ort hogonal to " }{TEXT 275 1 "v" }{XPPEDIT 18 0 "``[2]" "6#&%!G6#\"\"#" } {TEXT -1 5 " and " }{TEXT 276 1 "v" }{XPPEDIT 18 0 "``[3]" "6#&%!G6#\" \"$" }{TEXT -1 9 ". Thus, " }{TEXT 277 1 "v" }{XPPEDIT 18 0 "``[1]" " 6#&%!G6#\"\"\"" }{TEXT -1 41 " is not orthogonal to the null space, so " }{TEXT 278 1 "v" }{XPPEDIT 18 0 "``[1]" "6#&%!G6#\"\"\"" }{TEXT -1 4 " is " }{TEXT 288 3 "not" }{TEXT -1 21 " in the row space of " } {TEXT 418 1 "A" }{TEXT -1 14 ". The vector " }{TEXT 279 1 "v" } {XPPEDIT 18 0 "``[1]" "6#&%!G6#\"\"\"" }{TEXT -1 37 " has a component \+ in the row space of " }{TEXT 419 1 "A" }{TEXT -1 51 ", but it also has a component in the null space of " }{TEXT 420 1 "A" }{TEXT -1 1 "." } }{PARA 0 "" 0 "" {TEXT -1 0 "" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 39 "dotprod(v[1],v[2]);\ndotprod(v[1],v[3]);" }}}{EXCHG {PARA 0 "> \+ " 0 "" {MPLTEXT 1 0 0 "" }}}}{SECT 1 {PARA 3 "" 0 "" {TEXT -1 23 "Mini mum Length Solution" }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{EXCHG {PARA 0 "" 0 "" {TEXT -1 5 "When " }{TEXT 421 1 "A" }{TEXT -1 253 " has less t han maximal rank, the solution to the least-squares problem is not uni que. One way to select a unique solution from the resulting family of least-squares solutions is to pick the least-squares solution of mini mum length. Therefore, minimize" }}{PARA 0 "" 0 "" {TEXT -1 0 "" }} {PARA 267 "" 0 "" {TEXT -1 7 "|| X ||" }{XPPEDIT 18 0 "``^2" "6#*$%!G \"\"#" }{TEXT -1 3 " = " }{XPPEDIT 18 0 "1430*c[2]^2+1185*c[1]^2-2600* c[2]*c[1]+``(723500/13979)*c[2]-``(658488/13979)*c[1]+91593329/1954124 41" "6#,.*&\"%I9\"\"\"*$&%\"cG6#\"\"#F+F&F&*&\"%&=\"F&*$&F)6#F&F+F&F&* (\"%+EF&&F)6#F+F&&F)6#F&F&!\"\"*&-%!G6#*&\"'+NsF&\"&zR\"F7F&&F)6#F+F&F &*&-F:6#*&\"')[e'F&F>F7F&&F)6#F&F&F7*&\")HLf\"*F&\"*TCT&>F7F&" }{TEXT -1 1 " " }}{PARA 0 "" 0 "" {TEXT -1 8 "that is," }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 36 "q := simpli fy(norm(X,2),symbolic)^2;" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 35 "with respect to the two para meters " }{XPPEDIT 18 0 "c[1]" "6#&%\"cG6#\"\"\"" }{TEXT -1 5 " and " }{XPPEDIT 18 0 "c[2]" "6#&%\"cG6#\"\"#" }{TEXT -1 129 ". The usual te chniques of multivariable calculus, namely, equating partial derivativ es to zero and solving simultaneously, yield" }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 268 "" 0 "" {TEXT -1 1 " " }{XPPEDIT 18 0 "c[1] = 4184 /489265" "6#/&%\"cG6#\"\"\"*&\"%%=%F'\"'l#*[!\"\"" }{TEXT -1 4 " , " }{XPPEDIT 18 0 "c[2] = -13131/1272089" "6#/&%\"cG6#\"\"#,$*&\"&JJ\"\" \"\"\"(*3s7!\"\"F-" }{TEXT -1 8 " , and " }{TEXT 422 1 "X" }{XPPEDIT 18 0 "``^`+`" "6#)%!G%\"+G" }{TEXT -1 3 " = " }{XPPEDIT 18 0 "MATRIX([ [-13131/1272089], [4184/489265], [1411/6360445], [-939/1272089]])" "6# -%'MATRIXG6#7&7#,$*&\"&JJ\"\"\"\"\"(*3s7!\"\"F-7#*&\"%%=%F+\"'l#*[F-7# *&\"%69F+\"(X/O'F-7#,$*&\"$R*F+F,F-F-" }{TEXT -1 1 " " }}{PARA 0 "" 0 "" {TEXT -1 38 "that is, using the calculus, we obtain" }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 53 "qq := \+ solve(\{diff(q,c[1]),diff(q,c[2])\},\{c[1],c[2]\});" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 42 "a nd substituting the optimizing values of " }{XPPEDIT 18 0 "c[1]" "6#&% \"cG6#\"\"\"" }{TEXT -1 5 " and " }{XPPEDIT 18 0 "c[2]" "6#&%\"cG6#\" \"#" }{TEXT -1 2 " " }{TEXT 265 1 "X" }{TEXT -1 94 " to determine the least-squares solution of minimal length. Texts typically call this \+ vector " }{TEXT 266 1 "X" }{XPPEDIT 18 0 "``^`+`" "6#)%!G%\"+G" } {TEXT -1 30 " which, in Maple, we will call" }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 24 "Xplus := subs(qq, op(X));" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 27 "Adjusting the component of " }{TEXT 280 1 "X" }{TEXT -1 22 " in the null space of " }{TEXT 423 1 "A" }{TEXT -1 27 " to minimize the length of " }{TEXT 281 1 "X" }{TEXT -1 29 " is equivalent to projecting " }{TEXT 282 1 "X" }{TEXT -1 23 " onto the r ow space of " }{TEXT 424 1 "A" }{TEXT -1 24 ". But the only part of \+ " }{TEXT 283 1 "X" }{TEXT -1 40 " having a component in the row space \+ of " }{TEXT 425 1 "A" }{TEXT -1 4 " is " }{TEXT 284 1 "v" }{XPPEDIT 18 0 "``[1]" "6#&%!G6#\"\"\"" }{TEXT -1 103 ". Hence, the minimizatio n process we just went through is nothing more than finding the projec tion of " }{TEXT 285 1 "v" }{XPPEDIT 18 0 "``[1]" "6#&%!G6#\"\"\"" } {TEXT -1 23 " onto the row space of " }{TEXT 426 1 "A" }{TEXT -1 44 ". We can demonstrate this claim by forming " }{TEXT 427 1 "P" }{TEXT -1 50 "', a matrix which projects onto the row space of " }{TEXT 428 1 "A" }{TEXT -1 24 ". (Distinguish between " }{TEXT 429 1 "P" }{TEXT -1 6 "' and " }{TEXT 430 1 "P" }{TEXT -1 39 ", the matrix which projec ts vectors in " }{TEXT 431 1 "R" }{XPPEDIT 18 0 "``^r" "6#)%!G%\"rG" } {TEXT -1 3 " = " }{TEXT 432 1 "R" }{XPPEDIT 18 0 "``^5" "6#*$%!G\"\"& " }{TEXT -1 26 " onto the column space of " }{TEXT 433 1 "A" }{TEXT -1 33 ".) If the columns of the matrix " }{TEXT 434 1 "C" }{TEXT -1 34 " are a basis for the row space of " }{TEXT 435 1 "A" }{TEXT -1 6 " , then" }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 269 "" 0 "" {TEXT -1 1 " " }{XPPEDIT 18 0 "C=MATRIX([[0, 1], [1, 0], [-23/44, -5/11], [15/2 2, 7/11]])" "6#/%\"CG-%'MATRIXG6#7&7$\"\"!\"\"\"7$F+F*7$,$*&\"#BF+\"#W !\"\"F2,$*&\"\"&F+\"#6F2F27$*&\"#:F+\"#AF2*&\"\"(F+F6F2" }{TEXT -1 5 " , " }{TEXT 436 1 "P" }{TEXT -1 4 "' = " }{TEXT 437 1 "C" }{TEXT -1 1 "(" }{TEXT 438 1 "C" }{XPPEDIT 18 0 "``^T" "6#)%!G%\"TG" }{TEXT 439 1 "C" }{TEXT -1 1 ")" }{XPPEDIT 18 0 "``^`-1`" "6#)%!G%#-1G" }{TEXT 440 1 "C" }{XPPEDIT 18 0 "``^T" "6#)%!G%\"TG" }{TEXT -1 3 " = " } {XPPEDIT 18 0 "MATRIX([[673/910, -2/7, -17/91, 251/910], [-2/7, 24/35, -8/35, 2/7], [-17/91, -8/35, 93/455, -25/91], [251/910, 2/7, -25/91, \+ 337/910]])" "6#-%'MATRIXG6#7&7&*&\"$t'\"\"\"\"$5*!\"\",$*&\"\"#F*\"\"( F,F,,$*&\"# " 0 " " {MPLTEXT 1 0 1854 "g1 := plot([[-1,4],[1,4],[1,2],[-1,2],[-1,4]],col or=black):\ng2 := plot([[1,2],[2.5,2],[2.5,0],[1,0],[1,2]],color=black ):\ng3 := plot([[3,4],[4.5,4],[4.5,2],[3,2],[3,4]],color=black):\ng4 : = plot([[4.5,2],[6,2],[6,0],[4.5,0],[4.5,2]],color=black):\ng5 := arro w([1,3],[3,3],.05,.2,.1,color=black):\ng6 := arrow([2.5,.5],[4.5,2],.0 5,.2,.1,color=black):\ng7 := textplot(\{[0,4.2,`R`],[3.75,4.2,`R`]\},f ont=[TIMES,ROMAN,12]):\ng8 := textplot(\{[.15,4.4,`4`],[3.9,4.4,`5`]\} ,font=[TIMES,ROMAN,8]):\ng9 := textplot(\{[.7,3.8,`Row space of A`],[. 65,3.5,`dimension = 2`],[4,3.8,`Column`],[4.2,3.6,`space of A`],[4.2,3 .3,`dimension`],[4,3.05,`= 2 `],[2.25,.5,`Null space`],[1.9,.2,`of A` ],[2.5,-.2,`dimension = 2`]\},font=[TIMES,ROMAN,8], align=LEFT):\ng10 \+ := textplot(\{[4.8,.6,`Null space`],[5,.2,`of A`],[5.4,.3,`T`],[4.58,- .2,`dimension = 3`]\},font=[TIMES,ROMAN,8],align=RIGHT):\ng11 := arrow ([1,2],[-.6,2.6],.05,.2,.1,color=yellow):\ng12 := arrow([1,2],[1.3,1], .05,.2,.1,color=red):\ng13 := arrow([1,2],[1.7,1.4],.05,.2,.1,color=re d):\ng14 := arrow([1,2],[2,2.5],.05,.2,.1,color=green):\ng15 := arrow( [1,2],[1.6,3.7],.05,.2,.1,color=black):\ng16 := arrow([4.5,2],[3.3,2.6 ],.05,.2,.1,color=magenta):\ng17 := arrow([4.5,2],[4.9,3.7],.05,.2,.1, color=blue):\ng18 := arrow([4.5,2],[4.8,1],.05,.2,.1,color=orange):\ng 19 := textplot(\{[-.4,2.9,`X`],[1.85,3.6,`X`], [2.1,2.5,`v`], [1.8,1.6 ,`v`],[1.5,1,`v`],[5.2,3.6,`y`],[5,1.25,`y`],[3.9,2.7,`y`], [.5,2.9,`v `]\}, font=[TIMES,BOLD,10]):\ng20 := textplot(\{[2.2,2.4,`1`],[.6,2.8, `1`],[1.93,1.5,`2`],[1.63,.9,`3`], [-.18,3,`+`]\}, font=[TIMES,ROMAN,8 ]):\ng21 := textplot(\{[.15,2.9,`= P'`],[3.7,2.7,`P`]\},font=[TIMES,RO MAN,10]):\ng22 := textplot([5.1,1.1,`_|_`],font=[TIMES,ROMAN,5]):\ng23 := textplot([5.45,1.1,`(P )`],font=[TIMES,ROMAN,6]):\ng24 := textpl ot([5.5,1.1,`y`],font=[TIMES,BOLD,6]):\ndisplay([g||(1..24)],scaling=c onstrained, axes=none);" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 " " }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 31 "The basis for the row space o f " }{XPPEDIT 18 0 "A" "6#%\"AG" }{TEXT -1 25 " can be obtained Maple \+ as" }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 19 "rsA := rowspace(A);" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 209 "As in Sectio n 34.5, the order of the vectors is not fixed. To be sure our calcula tions here agree with what is written in the text, we note that the ve ctor with the larger 1-norm is on the left in the matrix " }{XPPEDIT 18 0 "C" "6#%\"CG" }{TEXT -1 29 " in the text. Hence, we for " }} {PARA 0 "" 0 "" {TEXT -1 0 "" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 66 "C := augment(op(sort([op(rsA)], (a,b)->is(norm(a,1)>norm(b,1)))) );" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" }}}{EXCHG {PARA 0 " " 0 "" {TEXT -1 5 "From " }{TEXT 446 1 "C" }{TEXT -1 7 ", form " } {TEXT 447 1 "P" }{TEXT -1 4 "' = " }{TEXT 448 1 "C" }{TEXT -1 1 "(" } {TEXT 449 1 "C" }{XPPEDIT 18 0 "``^T" "6#)%!G%\"TG" }{TEXT 450 1 "C" } {TEXT -1 1 ")" }{XPPEDIT 18 0 "``^`-1`" "6#)%!G%#-1G" }{TEXT 451 1 "C " }{XPPEDIT 18 0 "``^T" "6#)%!G%\"TG" }{TEXT -1 80 ", the matrix which projects vectors onto the subspace spanned by the columns of " } {TEXT 452 1 "C" }{TEXT -1 33 ", that is, onto the row space of " } {TEXT 453 1 "A" }{TEXT -1 1 "." }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}} {EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 63 "`P'` := evalm(C &* inverse(t ranspose(C) &* C) &* transpose(C));" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 11 "The matrix " }{TEXT 454 1 "P" }{TEXT -1 17 "' should project " }{TEXT 294 1 "v" } {XPPEDIT 18 0 "``[1]" "6#&%!G6#\"\"\"" }{TEXT -1 5 " and " }{TEXT 295 1 "X" }{TEXT -1 60 " (the general least-squares solution) onto the row space of " }{TEXT 455 1 "A" }{TEXT -1 29 ". This projection should b e " }{TEXT 286 1 "X" }{XPPEDIT 18 0 "``^`+`" "6#)%!G%\"+G" }{TEXT -1 68 ", the least-squares solution of shortest length. Thus, the vector s " }{TEXT 456 1 "P" }{TEXT -1 1 "'" }{TEXT 296 1 "v" }{XPPEDIT 18 0 " ``[1]" "6#&%!G6#\"\"\"" }{TEXT -1 2 ", " }{TEXT 457 1 "P" }{TEXT -1 1 "'" }{TEXT 297 1 "X" }{TEXT -1 6 ", and " }{TEXT 298 1 "X" }{XPPEDIT 18 0 "``^`+`" "6#)%!G%\"+G" }{TEXT -1 64 " should all be the same, as \+ the following calculuations confirm." }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 52 "print(evalm(`P'` &* v[1]) , evalm(`P'` &* X), Xplus);" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" }}}}{SECT 1 {PARA 3 "" 0 "" {TEXT -1 13 "Pseudoinverse" }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{EXCHG {PARA 0 "" 0 "" {TEXT -1 45 "The over -determined system we are solving is " }{TEXT 458 1 "A" }{TEXT 299 1 " x" }{TEXT -1 3 " = " }{TEXT 319 1 "y" }{TEXT -1 67 ". There is no sol ution - the equations are inconsistent - because " }{TEXT 320 1 "y" } {TEXT -1 31 " is not in the column space of " }{TEXT 459 1 "A" }{TEXT -1 43 ". We settled for a least-squares solution " }{TEXT 300 1 "X" } {TEXT -1 109 " which was not unique. To extract a unique solution fro m the two-parameter family of solutions contained in " }{TEXT 301 1 "X " }{TEXT -1 11 ", we found " }{TEXT 302 1 "X" }{XPPEDIT 18 0 "``^`+`" "6#)%!G%\"+G" }{TEXT -1 116 ", the least-squares solution of minimum l ength, and confirmed that this unique solution resides in the row spac e of " }{TEXT 460 1 "A" }{TEXT -1 27 ", and is the projection of " } {TEXT 303 1 "X" }{TEXT -1 21 " onto that row space." }}{PARA 0 "" 0 " " {TEXT -1 0 "" }}{PARA 0 "" 0 "" {TEXT -1 14 "If the system " }{TEXT 461 1 "A" }{TEXT 304 1 "x" }{TEXT -1 3 " = " }{TEXT 321 1 "y" }{TEXT -1 72 " were a square system of maximal rank, there would be an invers e matrix " }{TEXT 462 1 "A" }{XPPEDIT 18 0 "``^`-1`" "6#)%!G%#-1G" } {TEXT -1 51 ", multiplication by which would yield the solution " } {TEXT 307 1 "x" }{TEXT -1 3 " = " }{TEXT 463 1 "A" }{XPPEDIT 18 0 "``^ `-1`" "6#)%!G%#-1G" }{TEXT 322 1 "y" }{TEXT -1 5 " via " }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 258 "" 0 "" {TEXT -1 1 " " }{TEXT 464 1 "A " }{XPPEDIT 18 0 "``^`-1`" "6#)%!G%#-1G" }{TEXT 465 1 "A" }{TEXT -1 1 " " }{TEXT 305 1 "x" }{TEXT -1 3 " = " }{TEXT 466 1 "A" }{XPPEDIT 18 0 "``^`-1`" "6#)%!G%#-1G" }{TEXT 323 1 "y" }{TEXT -1 3 " = " }{TEXT 306 2 "x " }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 0 "" 0 "" {TEXT -1 47 "For the over-determined case, we seek a matrix " }{TEXT 467 1 "A" }{XPPEDIT 18 0 "``^`+`" "6#)%!G%\"+G" }{TEXT -1 12 " with which " } {TEXT 324 1 "y" }{TEXT -1 78 " could be multiplied so as to yield the \+ minimum-length least-squares solution " }{TEXT 308 1 "X" }{XPPEDIT 18 0 "``^`+`" "6#)%!G%\"+G" }{TEXT -1 3 " = " }{TEXT 468 1 "A" }{XPPEDIT 18 0 "``^`+`" "6#)%!G%\"+G" }{TEXT 325 1 "y" }{TEXT -1 34 ". Such a m atrix will be called a " }{TEXT 309 13 "pseudoinverse" }{TEXT -1 67 ", and we now give two prescriptions for constructing such a matrix." }} {PARA 0 "" 0 "" {TEXT -1 0 "" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" }}}{PARA 0 "" 0 "" {TEXT -1 0 "" }}}{SECT 1 {PARA 3 "" 0 "" {TEXT -1 32 "Pseudoinverse by Curious Formula" }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{EXCHG {PARA 0 "" 0 "" {TEXT -1 46 "We seek a matrix, \+ the application of which to " }{TEXT 326 1 "y" }{TEXT -1 49 ", yields \+ the minimum-norm least-squares solution " }{TEXT 267 1 "X" }{XPPEDIT 18 0 "``^`+`" "6#)%!G%\"+G" }{TEXT -1 14 ". The matrix " }{TEXT 469 1 "A" }{XPPEDIT 18 0 "``^`+`" "6#)%!G%\"+G" }{TEXT -1 31 " can be comp uted by the formula" }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 259 "" 0 "" {TEXT -1 1 " " }{XPPEDIT 18 0 "A^`+`=limit([(A^T*A+s^2*I)^`-1`*A^T] ,s=0)" "6#/)%\"AG%\"+G-%&limitG6$7#*&),&*&)F%%\"TG\"\"\"F%F1F1*&%\"sG \"\"#%\"IGF1F1%#-1GF1)F%F0F1/F3\"\"!" }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 0 "" 0 "" {TEXT -1 54 "discussed briefly in [1]. When applied to the matrix " }{TEXT 475 1 "A" }{TEXT -1 25 " in Table 34.4, it yie lds" }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 270 "" 0 "" {TEXT -1 1 " \+ " }{XPPEDIT 18 0 "A^`+`=MATRIX([[-14211/12720890, -1418/908635, 41512/ 6360445, -2514/908635, 50311/12720890], [-1908/489265, -186/69895, 48/ 97853, 104/69895, -62/489265], [32391/12720890, 3817/1817270, -4100/12 72089, 436/908635, -11013/6360445], [-42867/12720890, -2551/908635, 28 544/6360445, -678/908635, 30917/12720890]])" "6#/)%\"AG%\"+G-%'MATRIXG 6#7&7',$*&\"&6U\"\"\"\"\")!*3s7!\"\"F1,$*&\"%=9F/\"'N'3*F1F1*&\"&7:%F/ \"(X/O'F1,$*&\"%9DF/F5F1F1*&\"&6.&F/F0F17',$*&\"%3>F/\"'l#*[F1F1,$*&\" $'=F/\"&&*)pF1F1*&\"#[F/\"&`y*F1*&\"$/\"F/FFF1,$*&\"#iF/FBF1F17'*&\"& \"RKF/F0F1*&\"% " 0 "" {MPLTEXT 1 0 81 "Aplus := m ap(limit,evalm(inverse(transpose(A) &* A + s^2) &* transpose(A)), s=0) ;" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 25 "A calculation then shows " }{TEXT 268 1 "X" } {XPPEDIT 18 0 "``^`+`" "6#)%!G%\"+G" }{TEXT -1 4 " = A" }{XPPEDIT 18 0 "``^`+`" "6#)%!G%\"+G" }{TEXT -1 1 " " }{TEXT 327 1 "y" }{TEXT -1 10 ", that is," }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}}{EXCHG {PARA 0 "> \+ " 0 "" {MPLTEXT 1 0 32 "print(Xplus, evalm(Aplus &* y));" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 39 "Moreover, we can also show the rank of " }{XPPEDIT 18 0 "A^`+` " "6#)%\"AG%\"+G" }{TEXT -1 10 " is 2, via" }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 12 "rank(Aplus);" }}} {EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 31 "Subsequently, we will show how " }{XPPEDIT 18 0 "A^`+`" " 6#)%\"AG%\"+G" }{TEXT -1 15 " is related to " }{XPPEDIT 18 0 "U*Sigma* V^T" "6#*(%\"UG\"\"\"%&SigmaGF%)%\"VG%\"TGF%" }{TEXT -1 38 ", the sing ular value decomposition of " }{TEXT 476 1 "A" }{TEXT -1 1 "." }} {PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 0 "" 0 "" {TEXT -1 11 "The matri x " }{TEXT 477 1 "A" }{XPPEDIT 18 0 "``^`+`" "6#)%!G%\"+G" }{TEXT -1 46 " behaves as an inverse on the column space of " }{TEXT 478 1 "A" } {TEXT -1 33 ", acting as follows. Vectors in " }{TEXT 479 1 "R" } {XPPEDIT 18 0 "``^r" "6#)%!G%\"rG" }{TEXT -1 3 " = " }{TEXT 480 1 "R" }{XPPEDIT 18 0 "``^5" "6#*$%!G\"\"&" }{TEXT -1 40 " are projected onto the column space of " }{TEXT 481 1 "A" }{TEXT -1 59 ", and those proj ections are \"pulled back\" to pre-images in " }{TEXT 482 1 "R" } {XPPEDIT 18 0 "``^c" "6#)%!G%\"cG" }{TEXT -1 3 " = " }{TEXT 483 1 "R" }{XPPEDIT 18 0 "``^4" "6#*$%!G\"\"%" }{TEXT -1 61 ". These pre-images are then projected onto the row space of " }{TEXT 484 1 "A" }{TEXT -1 12 ". A vector " }{TEXT 328 1 "y" }{TEXT -1 4 " in " }{TEXT 485 1 "R" }{XPPEDIT 18 0 "``^r" "6#)%!G%\"rG" }{TEXT -1 3 " = " }{TEXT 486 1 "R" }{XPPEDIT 18 0 "``^5" "6#*$%!G\"\"&" }{TEXT -1 55 " is projected onto the two-dimensional column space of " }{TEXT 487 1 "A" }{TEXT -1 52 ", and that projection is inverted to a pre-image in " }{TEXT 488 1 "R" }{XPPEDIT 18 0 "``^c" "6#)%!G%\"cG" }{TEXT -1 3 " = " } {TEXT 489 1 "R" }{XPPEDIT 18 0 "``^4" "6#*$%!G\"\"%" }{TEXT -1 74 " wh ich is not unique. By projecting this pre-image onto the row space of " }{TEXT 490 1 "A" }{TEXT -1 47 ", the unique, shortest least-squares solution, " }{TEXT 310 1 "X" }{XPPEDIT 18 0 "``^`+`" "6#)%!G%\"+G" } {TEXT -1 4 " = A" }{XPPEDIT 18 0 "``^`+`" "6#)%!G%\"+G" }{TEXT -1 1 " \+ " }{TEXT 329 1 "y" }{TEXT -1 13 " is obtained." }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" }}}{PARA 0 "" 0 "" {TEXT -1 0 "" }}}{SECT 1 {PARA 3 "" 0 "" {TEXT -1 22 "Pseudo inverse from SVD" }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{EXCHG {PARA 0 "" 0 "" {TEXT -1 39 "If the singular value decomposition of " }{TEXT 491 1 "A" }{TEXT -1 13 " is given by " }{XPPEDIT 18 0 "A=U*Sigma*V^T" "6#/ %\"AG*(%\"UG\"\"\"%&SigmaGF')%\"VG%\"TGF'" }{TEXT -1 31 ", the pseudoi nverse is given by" }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 271 "" 0 " " {TEXT -1 1 " " }{XPPEDIT 18 0 "A^`+`=V*Sigma^`+`*U^T" "6#/)%\"AG%\"+ G*(%\"VG\"\"\")%&SigmaGF&F))%\"UG%\"TGF)" }{TEXT -1 1 " " }}{PARA 0 " " 0 "" {TEXT -1 6 "where " }{XPPEDIT 18 0 "Sigma^`+`" "6#)%&SigmaG%\"+ G" }{TEXT -1 37 " is the transpose of that version of " }{XPPEDIT 18 0 "Sigma" "6#%&SigmaG" }{TEXT -1 82 " in which eqch nonzero singular v alue has been replaced by its reciprocal. Since " }{XPPEDIT 18 0 "Sig ma" "6#%&SigmaG" }{TEXT -1 7 ", like " }{TEXT 492 1 "A" }{TEXT -1 17 " , has dimensions " }{XPPEDIT 18 0 "r" "6#%\"rG" }{TEXT -1 1 " " } {TEXT 493 1 "x" }{TEXT -1 1 " " }{XPPEDIT 18 0 "c" "6#%\"cG" }{TEXT -1 5 " = 5 " }{TEXT 494 1 "x" }{TEXT -1 9 " 4, then " }{XPPEDIT 18 0 " Sigma^`+`" "6#)%&SigmaG%\"+G" }{TEXT -1 16 " has dimensions " } {XPPEDIT 18 0 "c" "6#%\"cG" }{TEXT -1 1 " " }{TEXT 495 1 "x" }{TEXT -1 1 " " }{XPPEDIT 18 0 "r" "6#%\"rG" }{TEXT -1 5 " = 4 " }{TEXT 496 1 "x" }{TEXT -1 29 " 5. For our example, we have" }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 272 "" 0 "" {TEXT -1 1 " " }{XPPEDIT 18 0 "Sigma ^`+`=MATRIX([[1/sqrt(20950+30*sqrt(204983)), 0, 0, 0, 0], [0, 1/sqrt(2 0950-30*sqrt(204983)), 0, 0, 0], [0, 0, 0, 0, 0], [0, 0, 0, 0, 0]])" " 6#/)%&SigmaG%\"+G-%'MATRIXG6#7&7'*&\"\"\"F--%%sqrtG6#,&\"&]4#F-*&\"#IF --F/6#\"'$)\\?F-F-!\"\"\"\"!F9F9F97'F9*&F-F--F/6#,&F2F-*&F4F--F/6#F7F- F8F8F9F9F97'F9F9F9F9F97'F9F9F9F9F9" }{TEXT -1 1 " " }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 0 "" 0 "" {TEXT -1 29 "and a calculation shows t hat " }{XPPEDIT 18 0 "A^`+`=V*Sigma^`+`*U^T" "6#/)%\"AG%\"+G*(%\"VG\" \"\")%&SigmaGF&F))%\"UG%\"TGF)" }{TEXT -1 1 "." }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 0 "" 0 "" {TEXT -1 78 "To implement such a calcu lation in Maple, we must have access to the matrices " }{XPPEDIT 18 0 "U, V" "6$%\"UG%\"VG" }{TEXT -1 6 ", and " }{XPPEDIT 18 0 "Sigma" "6#% &SigmaG" }{TEXT -1 96 ". To this end, we have copied and pasted the e xpressions for these matrices. Thus, the matrix " }{XPPEDIT 18 0 "V" "6#%\"VG" }{TEXT -1 12 " is given by" }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 775 "V := matrix([[739/2*(327 7+23*204983^(1/2))/(52054667238495+80926940385*204983^(1/2))^(1/2),-73 9/2*(-3277+23*204983^(1/2))/(52054667238495-80926940385*204983^(1/2))^ (1/2),0,1/910*910^(1/2)*237^(1/2)], [-1478*(2833+5*204983^(1/2))/ (520 54667238495+80926940385*204983^(1/2))^(1/2), 1478*(-2833+5*204983^(1/2 ))/(52054667238495-80926940385*204983^(1/2))^(1/2), 1/1185*1185^(1/2), 2/1659*910^(1/2)*237^(1/2)], [1638363/(52054667238495+80926940385*204 983^(1/2))^(1/2), 1638363/(52054667238495-80926940385*204983^(1/2))^(1 /2), 28/1185*1185^(1/2), 17/21567*910^(1/2)*237^(1/2)], [739/2*(-5641+ 204983^(1/2))/(52054667238495+80926940385*204983^(1/2))^(1/2),-739/2*( 5641+204983^(1/2))/(52054667238495-80926940385*204983^(1/2))^(1/2), 4/ 237*1185^(1/2), -251/215670*910^(1/2)*237^(1/2)]]);" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 11 "t he matris " }{XPPEDIT 18 0 "U" "6#%\"UG" }{TEXT -1 12 " is given by" } }{PARA 0 "" 0 "" {TEXT -1 0 "" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 1493 "U := matrix([[232785*(2833+5*204983^(1/2))/(20950+30*204983^(1 /2))^(1/2)/(52054667238495+80926940385*204983^(1/2))^(1/2), -232785*(- 2833+5*204983^(1/2))/(20950-30*204983^(1/2))^(1/2)/(52054667238495-809 26940385*204983^(1/2))^(1/2), 8/7273*7273^(1/2), -693/1576163*1517^(1/ 2)*1039^(1/2), 468/3029449*1997^(1/2)*1517^(1/2)], [3695*(114083+181*2 04983^(1/2))/(20950+30*204983^(1/2))^(1/2)/(52054667238495+80926940385 *204983^(1/2))^(1/2), -3695*(-114083+181*204983^(1/2))/(20950-30*20498 3^(1/2))^(1/2)/(52054667238495-80926940385*204983^(1/2))^(1/2), 0, 1/1 517*1517^(1/2)*1039^(1/2), 244/3029449*1997^(1/2)*1517^(1/2)], [29560* (4903+29*204983^(1/2))/(20950+30*204983^(1/2))^(1/2)/(52054667238495+8 0926940385*204983^(1/2))^(1/2), -29560*(-4903+29*204983^(1/2))/(20950- 30*204983^(1/2))^(1/2)/(52054667238495-80926940385*204983^(1/2))^(1/2) , 45/7273*7273^(1/2), 128/1576163*1517^(1/2)*1039^(1/2), 544/3029449*1 997^(1/2)*1517^(1/2)], [-29560*(12367+32*204983^(1/2))/(20950+30*20498 3^(1/2))^(1/2)/(52054667238495+80926940385*204983^(1/2))^(1/2), 29560* (-12367+32*204983^(1/2))/(20950-30*204983^(1/2))^(1/2)/(52054667238495 -80926940385*204983^(1/2))^(1/2), 0, 0, 1/1997*1997^(1/2)*1517^(1/2)], [22170*(7391+30*204983^(1/2))/(20950+30*204983^(1/2))^(1/2)/(52054667 238495+80926940385*204983^(1/2))^(1/2), -22170*(-7391+30*204983^(1/2)) /(20950-30*204983^(1/2))^(1/2)/(52054667238495-80926940385*204983^(1/2 ))^(1/2), -72/7273*7273^(1/2), 3/1576163*1517^(1/2)*1039^(1/2), 392/30 29449*1997^(1/2)*1517^(1/2)]]);" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 15 "and the matrix " } {XPPEDIT 18 0 "Sigma" "6#%&SigmaG" }{TEXT -1 12 " is given by" }} {PARA 0 "" 0 "" {TEXT -1 0 "" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 144 "Sigma := matrix([[(20950+30*204983^(1/2))^(1/2), 0, 0, 0], [0, \+ (20950-30*204983^(1/2))^(1/2), 0, 0], [0, 0, 0, 0], [0, 0, 0, 0], [0, \+ 0, 0, 0]]);" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 64 "Just to be sure the factors are correctly reproduced, test that " }{TEXT 497 1 "A" }{TEXT -1 3 " = " }{TEXT 498 1 "U" }{TEXT -1 1 " " }{XPPEDIT 18 0 "Sigma" "6#%&SigmaG" }{TEXT -1 1 " " }{TEXT 499 1 "V" }{XPPEDIT 18 0 "``^T" "6#)%!G%\"TG" }{TEXT -1 1 "." }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 59 "print(A, map(radnormal,evalm(U &* Sigma &* transpose( V))));" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 12 "The factors " }{TEXT 500 1 "U" }{TEXT -1 5 " an d " }{TEXT 501 1 "V" }{TEXT -1 111 " are the important ones; we also n eed the singular values themselves, so we again obtain them exactly wi th the " }{TEXT 269 12 "singularvals" }{TEXT -1 9 " command." }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 22 "s v := singularvals(A);" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" } }}{EXCHG {PARA 0 "" 0 "" {TEXT -1 44 "For our example, the simplest wa y to create " }{XPPEDIT 18 0 "Sigma^`+`" "6#)%&SigmaG%\"+G" }{TEXT -1 13 " is directly." }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}}{EXCHG {PARA 0 " > " 0 "" {MPLTEXT 1 0 49 "SigmaPlus := extend(diag(1/sv[4],1/sv[5]),2, 3,0);" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 26 "Next, compute the product " }{TEXT 502 1 "V" } {TEXT -1 1 " " }{XPPEDIT 18 0 "Sigma^`+`" "6#)%&SigmaG%\"+G" }{TEXT -1 1 " " }{TEXT 503 1 "U" }{XPPEDIT 18 0 "``^T" "6#)%!G%\"TG" }{TEXT -1 31 " and verify that it reproduces " }{TEXT 504 1 "A" }{XPPEDIT 18 0 "``^`+`" "6#)%!G%\"+G" }{TEXT -1 1 "." }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 104 "map(rationalize,evalm (V &* SigmaPlus &* transpose(U)));\nlprint(`and from before, we had`); \nprint(Aplus);" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" }}} {EXCHG {PARA 0 "" 0 "" {TEXT -1 32 "As predicted, the pseudoinverse " }{XPPEDIT 18 0 "A^`+`" "6#)%\"AG%\"+G" }{TEXT -1 73 " can be obtained \+ from the factors in the singular value decomposition of " }{TEXT 505 1 "A" }{TEXT -1 1 "." }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" }}}}{SECT 1 {PARA 3 "" 0 "" {TEXT -1 33 "Numeric SVD and the Pseudoinverse" }}{PARA 0 "" 0 "" {TEXT -1 0 "" }} {EXCHG {PARA 0 "" 0 "" {TEXT -1 111 "In practice, the singular value d ecomposition is computed numerically, not exactly. We therefore consi der how " }{TEXT 506 1 "A" }{XPPEDIT 18 0 "``^`+`" "6#)%!G%\"+G" } {TEXT -1 37 " might be assembled from the factors " }{XPPEDIT 18 0 "U[ f]" "6#&%\"UG6#%\"fG" }{TEXT -1 5 " and " }{XPPEDIT 18 0 "V[f]" "6#&% \"VG6#%\"fG" }{TEXT -1 65 " which were obtained in floating-point form in Section 34.5. If " }{XPPEDIT 18 0 "Sigma[f]" "6#&%&SigmaG6#%\"fG " }{TEXT -1 36 ", a numerically computed version of " }{XPPEDIT 18 0 " Sigma" "6#%&SigmaG" }{TEXT -1 37 " is obtained from the singular value s" }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 26 "SV := evalf(Svd(A,Uf,Vf));" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 32 "the recipr ocals which appear in " }{XPPEDIT 18 0 "Sigma^`+`" "6#)%&SigmaG%\"+G" }{TEXT -1 37 " would include a term with magnitude " }{XPPEDIT 18 0 "1 0^13" "6#*$\"#5\"#8" }{TEXT -1 182 ". Clearly, the numerically determ ined singular values must be filtered through a threshold barrier to e liminate values which are probably zero. If this is done with a thres hold of " }{XPPEDIT 18 0 "10^(-10)" "6#)\"#5,$F$!\"\"" }{TEXT -1 4 " a s " }{XPPEDIT 18 0 "Sigma[f]^`+`" "6#)&%&SigmaG6#%\"fG%\"+G" }{TEXT -1 21 " is constructed from " }{XPPEDIT 18 0 "Sigma" "6#%&SigmaG" } {TEXT -1 6 ", then" }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 273 "" 0 " " {TEXT -1 1 " " }{XPPEDIT 18 0 "Sigma[f]^`+`=MATRIX([[.5381283824e-2, 0, 0, 0, 0], [0, .1165038708e-1, 0, 0, 0], [0, 0, 0, 0, 0], [0, 0, 0, 0, 0]])" "6#/)&%&SigmaG6#%\"fG%\"+G-%'MATRIXG6#7&7'$\"+CQG\"Q&!#7\"\" !F2F2F27'F2$\"+3(Q];\"!#6F2F2F27'F2F2F2F2F27'F2F2F2F2F2" }{TEXT -1 1 " " }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 0 "" 0 "" {TEXT -1 7 "and i f " }{XPPEDIT 18 0 "A[f]^`+`=V[f]*``(Sigma[f]^`+`)*U[f]^T" "6#/)&%\"AG 6#%\"fG%\"+G*(&%\"VG6#F(\"\"\"-%!G6#)&%&SigmaG6#F(F)F.)&%\"UG6#F(%\"TG F." }{TEXT -1 6 ", then" }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 274 " " 0 "" {TEXT -1 1 " " }{XPPEDIT 18 0 "A[f]^`+`=MATRIX([[-.1117138816e- 2, -.1560582632e-2, .6526587369e-2, -.2766787544e-2, .3954990570e-2], \+ [-.3899727142e-2, -.2661134558e-2, .490531716e-3, .1487946205e-2, -.12 6720694e-3], [.2546284104e-2, .2100403352e-2, -.3223044929e-2, .479840 6398e-3, -.1731482624e-2], [-.3369811389e-2, -.2807507965e-2, .4487736 313e-2, -.7461742057e-3, .2430411708e-2]])" "6#/)&%\"AG6#%\"fG%\"+G-%' MATRIXG6#7&7',$$\"+;)Qr6\"!#7!\"\",$$\"+KEeg:F2F3$\"+pteElF2,$$\"+Wvym FF2F3$\"+q0*\\&RF27',$$\"+Urs**QF2F3,$$\"+eX8hEF2F3$\"*;<`!\\F2$\"+0i% z[\"F2,$$\"*%p?n7F2F37'$\"+/TGYDF2$\"+_LS+@F2,$$\"+H\\/BKF2F3$\"+)R1%) z%!#8,$$\"+CE[J " 0 "" {MPLTEXT 1 0 59 "Sigmaf := m atrix(5,4,(i,j) -> if i=j then SV[i] else 0 fi);" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 25 "Setti ng the threshold as " }{XPPEDIT 18 0 "10^`-10`" "6#)\"#5%$-10G" } {TEXT -1 11 ", we build " }{XPPEDIT 18 0 "Sigma[f]^`+`" "6#)&%&SigmaG6 #%\"fG%\"+G" }{TEXT -1 12 " directly as" }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 108 "SigmaPlusf := transpo se(matrix(5,4, (i,j) -> if i=j and Sigmaf[i,i]>10^(-10) then 1/Sigmaf[ i,i] else 0 fi));" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" }}} {EXCHG {PARA 0 "" 0 "" {TEXT -1 14 "We next build " }{XPPEDIT 18 0 "A[ f]^`+`=V[f]*``(Sigma[f]^`+`)*U[f]^T" "6#/)&%\"AG6#%\"fG%\"+G*(&%\"VG6# F(\"\"\"-%!G6#)&%&SigmaG6#F(F)F.)&%\"UG6#F(%\"TGF." }{TEXT -1 4 " as " }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 51 "Aplusf := evalm(Vf &* SigmaPlusf &* transpose(Uf));" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" }}}{EXCHG {PARA 0 "" 0 " " {TEXT -1 30 "Floating the exact version of " }{TEXT 511 1 "A" } {XPPEDIT 18 0 "``^`+`" "6#)%!G%\"+G" }{TEXT -1 8 " we find" }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 17 "eva lf(op(Aplus));" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" }}} {EXCHG {PARA 0 "" 0 "" {TEXT -1 20 "Finally, we compute " }{TEXT 513 1 "X" }{XPPEDIT 18 0 "``[f]^`+`" "6#)&%!G6#%\"fG%\"+G" }{TEXT -1 16 " \+ by converting " }{TEXT 514 1 "X" }{XPPEDIT 18 0 "``^`+`" "6#)%!G%\"+G " }{TEXT -1 46 " to floating-point form, and we compare it to " } {XPPEDIT 18 0 "``(A[f]^`+`)" "6#-%!G6#)&%\"AG6#%\"fG%\"+G" }{TEXT -1 1 " " }{TEXT 512 1 "y" }{XPPEDIT 18 0 "``[f]" "6#&%!G6#%\"fG" }{TEXT -1 14 " , as follows." }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 54 "print(evalf(op(Xplus)),evalm(Aplusf &* evalf(op(y))));" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" }}} }{SECT 1 {PARA 3 "" 0 "" {TEXT -1 10 "References" }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 0 "" 0 "" {TEXT -1 96 "[1] Richard Pavelle, edit or, Applications of Computer Algebra, Kluwer Academic Publishers, 1985 ." }}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" }}}}}{MARK "1" 0 } {VIEWOPTS 1 1 0 1 1 1803 1 1 1 1 }{PAGENUMBERS 0 1 2 33 1 1 }