From help-octave-request at bevo dot che dot wisc dot edu Thu Oct 26 14:26:39 2000 Subject: Iterative Proportional Fitting From: Stefan Hrafn Jonsson To: help-octave at bevo dot che dot wisc dot edu Date: Thu, 26 Oct 2000 15:26:33 -0400 (EDT) Hi all I am new to Octave I am trying to find a way to use Iterative Proportional Fitting to solve simoltanioysly 13 equations. I do have a Maple code to do this but I would really like to do this in Octave, partly because I have not been able to find Maple for unix on this (Penn State) campus (except too old versions) and bacause people around tell me Octave should be good for this task. So far I have been able to construct in octave all the given values needed to solve the equations. I am stuck in the octave code wehre I go about to solve equations. To demostrate what equation I want I provide the Maple code used to solve this. Msd = 0.00068 Msm = 0.13446 Mmd = 0.00037 Mmw = 0.00116 Mmv = 0.06563 Mwd = 0.00148 Mwm = 0.03334 Mvd = 0.00113 Mvm = 0.5128 Ms = [Msd+Msm,0,0,0;-Msm,Mmd+Mmw+Mmv,-Mwm,-Mvm; 0,-Mmw,Mwd+Mwm,0; 0,-Mmv,0,Mvd+Mvm ] II = [1,0,0,0;0,1,0,0;0,0,1,0;0,0,0,1] Msplusi = II + 0.5*Ms Msplusinv = inverse(Msplusi) Msminusi = II - 0.5*Ms Spi = Msminusi * Msplusinv # Observed rates for from 1995 US MSLT, Females 22-23 aMsd=0.00058 aMsm=0.09772 aMmd=0.0003 aMmw=0.00112 aMmv=0.04588 aMwd=0.0014 aMwm=0.05119 aMvd=0.000753 aMvm=0.22632 AM = [aMsd+aMsm,0,0,0;-aMsm,aMmd+aMmw+aMmv,-aMwm,-aMvm;0,-aMmw,aMwd+aMwm,0;0,-aMmv,0,aMvd+aMvm] AMplusi = II + 0.5*AM AMplusinv =inverse(AMplusi) AMminusi = II - 0.5*AM Api = AMminusi*AMplusinv x = [76411;20716;64;1814] y = [69250;27035;85;2432] Cpi = Api * x Sctot1 = Spi(1,1)+Spi(2,1)+Spi(3,1)+Spi(4,1) Sctot2 = Spi(1,2)+Spi(2,2)+Spi(3,2)+Spi(4,2) Sctot3 = Spi(1,3)+Spi(2,3)+Spi(3,3)+Spi(4,3) Sctot4 = Spi(1,4)+Spi(2,4)+Spi(3,4)+Spi(4,4) # Now how do we solve the 13 equations #function Z = f (x) # Z = [z11,0,0,0;z21,z22,z23,z24;z31,z32,z33,z34;z41,z42,z43,z44] #endfunction # #The Maple way of diong this. Z := matrix(4,4,[z11,0,0,0i],[z21,z22,z23,z24],[z31,z32,z33,z34],[z41,z42,z43,z44]]); solve ({ y[1,1] = z11*x[1,1], y[2,1] = z21*x[1,1]+z22*z[2,1]+z23*x[3,1]+z24*x[4,1], y[3,1] = z31*x[1,1]+z32*z[2,1]+z33*x[3,1]+z34*x[4,1], y[4,1] = z41*x[1,1]+z42*z[2,1]+z43*x[3,1]+z44*x[4,1], Sctot2/Sctot1 = (z32+z32+z42)/(11+z21+z31+z41), Sctot3/Sctot1 = (z33+z33+z43)/(11+z21+z31+z41), Sctot4/Sctot1 = (z34+z34+z44)/(11+z21+z31+z41), z21*z42/(z22*z41) = Spi[2,1]*Spi[4,2]/(Spi[2,2]*Spi[4,1]), z21*z43/(z23*z41) = Spi[2,1]*Spi[4,3]/(Spi[2,3]*Spi[4,1]), z21*z44/(z24*z41) = Spi[2,1]*Spi[4,4]/(Spi[2,4]*Spi[4,1]), z31*z42/(z32*z41) = Spi[3,1]*Spi[4,2]/(Spi[2,2]*Spi[4,1]), z31*z43/(z33*z41) = Spi[3,1]*Spi[4,3]/(Spi[2,3]*Spi[4,1]), z31*z44/(z34*z41) = Spi[3,1]*Spi[4,4]/(Spi[2,4]*Spi[4,1])}, {z11,z21,z31,z41,z22,z32,z42,z23,z33,z43,z24,z34,z44}) Any help or hints will be apreciated. Thanks Stefan _________________________ Stefan Hrafn Jonsson 701 Oswald Tower Department of Sociology and The Population Research Institute The Pennsylvania State University, University Park, PA 16802 _________________________ ------------------------------------------------------------- Octave is freely available under the terms of the GNU GPL. Octave's home on the web: http://www.octave.org How to fund new projects: http://www.octave.org/funding.html Subscription information: http://www.octave.org/archive.html -------------------------------------------------------------